CMIP6 prep

---

library(raster)

library(plyr)
library(dplyr)

library(ggplot2)
library(viridisLite)
library(colorspace)
library(RColorBrewer)

library(sf)
library(stringr)

library(ggpubr)
library(kableExtra)

memory.limit(size = 160000)
## [1] Inf

#display.brewer.pal(n = 11, name ="RdYlGn")
colorvec <- brewer.pal(n = 11, name ="RdYlGn")
col.cat <- c("Hyperarid" = colorvec[2], "Arid"= colorvec[4],"Semi-arid"= colorvec[5], "Dry subhumid" = colorvec[8], "Humid"= colorvec[10], "Cold" = "powderblue")

1 Data

land_mask <- raster("Masks/land_sea_mask_1degree.nc4") 
land_mask.df <- as.data.frame(land_mask, xy = T) %>% setNames(c("lon","lat","lm")) # 0 if ocean, 1 if land


elev <- raster("Worldclim/wc2.1_10m/wc2.1_10m_elev.tif") 
elev.df <- elev %>% projectRaster(to = land_mask) %>% as.data.frame(xy = T) %>% setNames(c("lon","lat", "z"))


ipcc_regions <- shapefile("Masks/IPCC-WGI-reference-regions-v4.shp") %>% spTransform(crs("EPSG:4326")) 
ipcc_regions.raster <- ipcc_regions %>% rasterize(land_mask)
ipcc_regions.df <- as.data.frame(ipcc_regions.raster, xy = T) %>% setNames(c("lon","lat","Continent","Type","Name","Acronym"))

tas.all_annual <- read.table("CMIP6/tas/tas.all_annual.txt")
pr.all_annual <- read.table("CMIP6/pr/pr.all_annual.txt")
sfcWind.all_annual <- read.table("CMIP6/sfcWind/sfcWind.all_annual.txt")
hfls.all_annual <- read.table("CMIP6/hfls/hfls.all_annual.txt")
hfss.all_annual <- read.table("CMIP6/hfss/hfss.all_annual.txt")
hurs.all_annual <- read.table("CMIP6/hurs/hurs.all_annual.txt")

cmip6_annual <- tas.all_annual %>% merge(pr.all_annual, by = c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_annual, by = c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_annual, by = c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_annual, by = c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_annual, by = c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_annual, "cmip6_annual.txt")

tas.all_january <- read.table("CMIP6/tas/tas.all_january.txt")
pr.all_january <- read.table("CMIP6/pr/pr.all_january.txt")
sfcWind.all_january <- read.table("CMIP6/sfcWind/sfcWind.all_january.txt")
hfls.all_january <- read.table("CMIP6/hfls/hfls.all_january.txt")
hfss.all_january <- read.table("CMIP6/hfss/hfss.all_january.txt")
hurs.all_january <- read.table("CMIP6/hurs/hurs.all_january.txt")

cmip6_january <- tas.all_january %>% merge(pr.all_january, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_january, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_january, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_january, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_january, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_january, "cmip6_january.txt")

tas.all_february <- read.table("CMIP6/tas/tas.all_february.txt")
pr.all_february <- read.table("CMIP6/pr/pr.all_february.txt")
sfcWind.all_february <- read.table("CMIP6/sfcWind/sfcWind.all_february.txt")
hfls.all_february <- read.table("CMIP6/hfls/hfls.all_february.txt")
hfss.all_february <- read.table("CMIP6/hfss/hfss.all_february.txt")
hurs.all_february <- read.table("CMIP6/hurs/hurs.all_february.txt")

cmip6_february <- tas.all_february %>% merge(pr.all_february, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_february, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_february, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_february, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_february, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_february, "cmip6_february.txt")

tas.all_march <- read.table("CMIP6/tas/tas.all_march.txt")
pr.all_march <- read.table("CMIP6/pr/pr.all_march.txt")
sfcWind.all_march <- read.table("CMIP6/sfcWind/sfcWind.all_march.txt")
hfls.all_march <- read.table("CMIP6/hfls/hfls.all_march.txt")
hfss.all_march <- read.table("CMIP6/hfss/hfss.all_march.txt")
hurs.all_march <- read.table("CMIP6/hurs/hurs.all_march.txt")

cmip6_march <- tas.all_march %>% merge(pr.all_march, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_march, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_march, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_march, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_march, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_march, "cmip6_march.txt")

tas.all_april <- read.table("CMIP6/tas/tas.all_april.txt")
pr.all_april <- read.table("CMIP6/pr/pr.all_april.txt")
sfcWind.all_april <- read.table("CMIP6/sfcWind/sfcWind.all_april.txt")
hfls.all_april <- read.table("CMIP6/hfls/hfls.all_april.txt")
hfss.all_april <- read.table("CMIP6/hfss/hfss.all_april.txt")
hurs.all_april <- read.table("CMIP6/hurs/hurs.all_april.txt")

cmip6_april <- tas.all_april %>% merge(pr.all_april, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_april, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_april, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_april, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_april, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_april, "cmip6_april.txt")

tas.all_may <- read.table("CMIP6/tas/tas.all_may.txt")
pr.all_may <- read.table("CMIP6/pr/pr.all_may.txt")
sfcWind.all_may <- read.table("CMIP6/sfcWind/sfcWind.all_may.txt")
hfls.all_may <- read.table("CMIP6/hfls/hfls.all_may.txt")
hfss.all_may <- read.table("CMIP6/hfss/hfss.all_may.txt")
hurs.all_may <- read.table("CMIP6/hurs/hurs.all_may.txt")

cmip6_may <- tas.all_may %>% merge(pr.all_may, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_may, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_may, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_may, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_may, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_may, "cmip6_may.txt")

tas.all_june <- read.table("CMIP6/tas/tas.all_june.txt")
pr.all_june <- read.table("CMIP6/pr/pr.all_june.txt")
sfcWind.all_june <- read.table("CMIP6/sfcWind/sfcWind.all_june.txt")
hfls.all_june <- read.table("CMIP6/hfls/hfls.all_june.txt")
hfss.all_june <- read.table("CMIP6/hfss/hfss.all_june.txt")
hurs.all_june <- read.table("CMIP6/hurs/hurs.all_june.txt")

cmip6_june <- tas.all_june %>% merge(pr.all_june, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_june, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_june, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_june, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_june, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_june, "cmip6_june.txt")

tas.all_july <- read.table("CMIP6/tas/tas.all_july.txt")
pr.all_july <- read.table("CMIP6/pr/pr.all_july.txt")
sfcWind.all_july <- read.table("CMIP6/sfcWind/sfcWind.all_july.txt")
hfls.all_july <- read.table("CMIP6/hfls/hfls.all_july.txt")
hfss.all_july <- read.table("CMIP6/hfss/hfss.all_july.txt")
hurs.all_july <- read.table("CMIP6/hurs/hurs.all_july.txt")

cmip6_july <- tas.all_july %>% merge(pr.all_july, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_july, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_july, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_july, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_july, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_july, "cmip6_july.txt")

tas.all_august <- read.table("CMIP6/tas/tas.all_august.txt")
pr.all_august <- read.table("CMIP6/pr/pr.all_august.txt")
sfcWind.all_august <- read.table("CMIP6/sfcWind/sfcWind.all_august.txt")
hfls.all_august <- read.table("CMIP6/hfls/hfls.all_august.txt")
hfss.all_august <- read.table("CMIP6/hfss/hfss.all_august.txt")
hurs.all_august <- read.table("CMIP6/hurs/hurs.all_august.txt")

cmip6_august <- tas.all_august %>% merge(pr.all_august, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_august, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_august, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_august, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_august, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_august, "cmip6_august.txt")

tas.all_september <- read.table("CMIP6/tas/tas.all_september.txt")
pr.all_september <- read.table("CMIP6/pr/pr.all_september.txt")
sfcWind.all_september <- read.table("CMIP6/sfcWind/sfcWind.all_september.txt")
hfls.all_september <- read.table("CMIP6/hfls/hfls.all_september.txt")
hfss.all_september <- read.table("CMIP6/hfss/hfss.all_september.txt")
hurs.all_september <- read.table("CMIP6/hurs/hurs.all_september.txt")

cmip6_september <- tas.all_september %>% merge(pr.all_september, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_september, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_september, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_september, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_september, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_september, "cmip6_september.txt")

tas.all_october <- read.table("CMIP6/tas/tas.all_october.txt")
pr.all_october <- read.table("CMIP6/pr/pr.all_october.txt")
sfcWind.all_october <- read.table("CMIP6/sfcWind/sfcWind.all_october.txt")
hfls.all_october <- read.table("CMIP6/hfls/hfls.all_october.txt")
hfss.all_october <- read.table("CMIP6/hfss/hfss.all_october.txt")
hurs.all_october <- read.table("CMIP6/hurs/hurs.all_october.txt")

cmip6_october <- tas.all_october %>% merge(pr.all_october, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_october, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_october, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_october, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_october, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_october, "cmip6_october.txt")

tas.all_november <- read.table("CMIP6/tas/tas.all_november.txt")
pr.all_november <- read.table("CMIP6/pr/pr.all_november.txt")
sfcWind.all_november <- read.table("CMIP6/sfcWind/sfcWind.all_november.txt")
hfls.all_november <- read.table("CMIP6/hfls/hfls.all_november.txt")
hfss.all_november <- read.table("CMIP6/hfss/hfss.all_november.txt")
hurs.all_november <- read.table("CMIP6/hurs/hurs.all_november.txt")

cmip6_november <- tas.all_november %>% merge(pr.all_november, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_november, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_november, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_november, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_november, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_november, "cmip6_november.txt")

tas.all_december <- read.table("CMIP6/tas/tas.all_december.txt")
pr.all_december <- read.table("CMIP6/pr/pr.all_december.txt")
sfcWind.all_december <- read.table("CMIP6/sfcWind/sfcWind.all_december.txt")
hfls.all_december <- read.table("CMIP6/hfls/hfls.all_december.txt")
hfss.all_december <- read.table("CMIP6/hfss/hfss.all_december.txt")
hurs.all_december <- read.table("CMIP6/hurs/hurs.all_december.txt")

cmip6_december <- tas.all_december %>% merge(pr.all_december, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(sfcWind.all_december, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfls.all_december, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hfss.all_december, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>%
  merge(hurs.all_december, by = c("lon","lat","model","period","month","lm", "Continent","Type","Name","Acronym","source")) %>% 
  merge(elev.df, by = c("lon", "lat"))

write.table(cmip6_december, "cmip6_december.txt")

2 Aridity Index

2.1 De Martonne

cmip6 <- read.table("cmip6_annual.txt")

breaks.martonne <- c(0,5,10,20,"30","40",Inf)
cat.martonne <- c("[0,5]" = "Desert", "(5,10]"= "Arid", "(10,20]" = "Semi-arid","(20,30]" = "Temperate", "(30,40]" = "Humid", "(40,Inf]" = "Forest")
col.martonne <- c("Desert" = colorvec[2], "Arid"= colorvec[4],"Semi-arid"= colorvec[5], "Temperate" = colorvec[8], "Humid"= colorvec[10],"Forest"= colorvec[11])

cmip6$AIm <- with(cmip6,pr*60*60*24*365/(tas-273.15+10)) # pr in mm/y, tas in Ceslsius


cmip6$cat.AIm <- cmip6$AIm %>% cut(breaks.martonne, include.lowest = T) %>% revalue(cat.martonne)


map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "1970_2000" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AIm))+
    scale_fill_manual(values = col.martonne, na.translate = F)+
    labs(title = paste(i, "1970-2000", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

2.2 UNEP

ea from Rh:

ea = 6.108(RH/100)exp(17.27*T/(T+237.3)) Refs https://www.nature.com/articles/s41598-023-40499-6 : Allen and Tetens

id_vars <- c("lon","lat","model","period","lm", "Continent","Type","Name","Acronym","source")

cmip6_annual <- mutate(read.table("cmip6_annual.txt"),
                       t = tas - 273.15,
                       Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                       P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                       es = 0.61078*exp(17.27*t/(t+237.3)),
                       ea = (hurs/100)*es,
                       delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                       gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                       ET0_annual = ((0.408*delta*Rn+gamma*900/((t)+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind)))) 

cmip6_01 <- mutate(read.table("cmip6_january.txt"),
                   t = tas - 273.15,
                   pr_01 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_01 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind)))) # Evapotranspiration, G considered to be negligible. in mm/day. Wind speed at 10 m is converted to wind speed at 2 m by multiplying by 0.748

cmip6_02 <- mutate(read.table("cmip6_february.txt"),
                   t = tas - 273.15,
                   pr_02 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_02 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_03 <- mutate(read.table("cmip6_march.txt"),
                   t = tas - 273.15,
                   pr_03 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_03 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_04 <- mutate(read.table("cmip6_april.txt"),
                   t = tas - 273.15,
                   pr_04 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_04 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_05 <- mutate(read.table("cmip6_may.txt"),
                   t = tas - 273.15,
                   pr_05 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_05 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_06 <- mutate(read.table("cmip6_june.txt"),
                   t = tas - 273.15,
                   pr_06 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_06 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_07 <- mutate(read.table("cmip6_july.txt"),
                   t = tas - 273.15,
                   pr_07 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_07 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_08 <- mutate(read.table("cmip6_august.txt"),
                   t = tas - 273.15,
                   pr_08 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_08 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_09 <- mutate(read.table("cmip6_september.txt"),
                   t = tas - 273.15,
                   pr_09 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_09 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_10 <- mutate(read.table("cmip6_october.txt"),
                   t = tas - 273.15,
                   pr_10 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_10 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_11 <- mutate(read.table("cmip6_november.txt"),
                   t = tas - 273.15,
                   pr_11 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_11 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

cmip6_12 <- mutate(read.table("cmip6_december.txt"),
                   t = tas - 273.15,
                   pr_12 = pr*60*60*24, # mm/day
                   Rn = (hfls + hfss) * 86400 * 1e-6, # conversion from W/m2 to MJ/m2/day.
                   P = 101.3 * ((293-0.0065*z)/293)^(5.26),
                   es = 0.61078*exp(17.27*t/(t+237.3)),
                   ea = (hurs/100)*es,
                   delta = 4098*es/(t+237.3)^2, #slope vapor pressure
                   gamma = 0.00163*P/(2.501-52.6361e-3*t), #psychrometric constant
                   ET0_12 = ((0.408*delta*Rn+gamma*900/(t+273)*0.748*sfcWind*(es-ea))/(delta+gamma*(1+0.34*0.748*sfcWind))))

et0_monthly <- select(cmip6_annual, id_vars) %>%
  merge(select(cmip6_01, c(id_vars, "ET0_01", "pr_01")), by = id_vars) %>%
  merge(select(cmip6_02, c(id_vars, "ET0_02", "pr_02")), by = id_vars) %>%
  merge(select(cmip6_03, c(id_vars, "ET0_03", "pr_03")), by = id_vars) %>%
  merge(select(cmip6_04, c(id_vars, "ET0_04", "pr_04")), by = id_vars) %>%
  merge(select(cmip6_05, c(id_vars, "ET0_05", "pr_05")), by = id_vars) %>%
  merge(select(cmip6_06, c(id_vars, "ET0_06", "pr_06")), by = id_vars) %>%
  merge(select(cmip6_07, c(id_vars, "ET0_07", "pr_07")), by = id_vars) %>%
  merge(select(cmip6_08, c(id_vars, "ET0_08", "pr_08")), by = id_vars) %>%
  merge(select(cmip6_09, c(id_vars, "ET0_09", "pr_09")), by = id_vars) %>%
  merge(select(cmip6_10, c(id_vars, "ET0_10", "pr_10")), by = id_vars) %>%   
  merge(select(cmip6_11, c(id_vars, "ET0_11", "pr_11")), by = id_vars) %>%
  merge(select(cmip6_12, c(id_vars, "ET0_12", "pr_12")), by = id_vars)

et0_monthly$spr <- with(et0_monthly, pr_01*31 + pr_02 * 28 + pr_03 * 31 + pr_04 * 30 + pr_05 * 31 + pr_06 * 30 +
                          pr_07 * 31 + pr_08 * 31 + pr_09 * 30 + pr_10 * 31 + pr_11 * 30 + pr_12 * 31)
et0_monthly$sET0 <- with(et0_monthly, ET0_01*31 + ET0_02 * 28 + ET0_03 * 31 + ET0_04 * 30 + ET0_05 * 31 + ET0_06 * 30 +
                           ET0_07 * 31 + ET0_08 * 31 + ET0_09 * 30 + ET0_10 * 31 + ET0_11 * 30 + ET0_12 * 31)

cmip6 <- merge(cmip6_annual, et0_monthly, by = id_vars)
cmip6$AI <- with(cmip6, abs(spr/sET0))  
cmip6$AI_annual <- with(cmip6, abs((pr*60*60*24*365)/(ET0_annual*365)))

breaks.unesco <- c(-Inf,0,0.03,0.2,0.5,0.65,Inf)
cat.unesco <- c("(0,0.03]"= "Hyperarid", "(0.03,0.2]" = "Arid", "(0.2,0.5]" = "Semi-arid", "(0.5,0.65]" = "Dry subhumid","(0.65, Inf]" = "Humid", "(-Inf,0]" = "Cold")
col.cat <- c("Hyperarid" = colorvec[2], "Arid"= colorvec[4],"Semi-arid"= colorvec[5], "Dry subhumid" = colorvec[8], "Humid"= colorvec[10], "Cold" = "powderblue")

cmip6$cat.AI <- cmip6$AI %>% cut(breaks = breaks.unesco) %>% revalue(cat.unesco)
cmip6$cat.AI[which(cmip6$sET0 < 400)] <- "Cold"
cmip6$cat.AI %>% unique()

cmip6$cat.AI <- cmip6$AI %>% cut(breaks = breaks.unesco) %>% 
  revalue(cat.unesco)


cmip6$cat.AI_annual <- cmip6$AI_annual %>% cut(breaks = breaks.unesco) %>% 
  revalue(cat.unesco)

cmip6$cat.AI_annual[which(cmip6$ET0_annual*365 < 400)] <- "Cold"

write.table(cmip6, "cmip6_AI.txt")

3 Maps AI by model, all % calculated without Antarctica

cmip6 <- read.table("cmip6_AI.txt")

breaks.unesco <- c(-Inf,0.03,0.2,0.5,0.65,Inf)
cat.unesco <- c("(-Inf,0.03]"= "Hyperarid", "(0.03,0.2]" = "Arid", "(0.2,0.5]" = "Semi-arid", "(0.5,0.65]" = "Dry subhumid","(0.65, Inf]" = "Humid")

3.1 1850-1880

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "1850_1880" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "1850-1880", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "1850_1880"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))


tab <- cmip6 %>% subset(period == "1850_1880"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)


knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	5.1	9.0	13.2	6.1	33.4	38.0	28.3	0.3
CESM	5.2	10.8	14.7	5.0	35.7	35.6	28.5	0.3
CMCC	5.6	7.2	9.5	5.3	27.6	42.1	30.1	0.3
CMCC-ESM2	6.0	8.2	7.5	2.8	24.5	24.0	51.2	0.3
CNRM	5.7	12.1	11.4	6.6	35.8	38.0	25.9	0.2
EC-Earth3	8.6	9.0	11.0	4.5	33.1	36.7	30.1	0.2
FGOALS	5.6	12.8	15.9	7.2	41.5	29.9	28.4	0.2
GFDL-ESM4	5.6	9.3	8.7	4.1	27.7	37.0	34.9	0.2
INM	4.4	8.2	13.4	5.8	31.8	40.6	27.4	0.3
INM-CM5	3.5	8.8	11.5	5.6	29.4	43.2	27.1	0.3
MPI	8.5	10.9	10.1	4.4	33.9	34.5	31.4	0.2
MRI	6.3	10.6	8.7	3.8	29.4	39.6	30.9	0.3
NorESM-2-MM	5.1	10.4	18.4	5.9	39.8	35.4	24.4	0.3
mean	5.8	9.8	11.8	5.2	32.6	36.5	30.7	0.3
sd	1.4	1.6	3.2	1.2	7.4	5.1	6.7	0.1

3.2 1970-2000

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "1970_2000" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "1850-1880", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

library(ggrepel) 

pie <- cmip6 %>% subset(period == "1970_2000" & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  select(c("source", "cat.AI")) %>%
  reshape2::melt(id.vars = c("source")) %>%
  group_by(source, value) %>%
  summarise(count = n()) %>%
  mutate(perc = round(count/sum(count)*100, 1),
         lab.y = rev(cumsum(rev(count))) - count*0.5) %>% 
  ungroup() 

pie_list <- list()

for(i in unique(pie$source)){
  p <- ggplot(subset(pie, source == i))+
    geom_bar(aes(x = "", y = count, fill = value), width = 1, stat = "identity")+
    coord_polar("y", start = 0)+
    scale_fill_manual(values = col.cat)+
    geom_label_repel(aes(x = "", y = lab.y, label = perc), nudge_x = 0.5)+
    labs(fill = "", title = i)+
    theme_void()+theme(legend.position = "right")
  pie_list[[i]] <- p
}

ggpubr::ggarrange(plotlist = pie_list, nrow = 4, ncol = 4, common.legend = T, legend = "bottom")

pie <- cmip6 %>% subset(period == "1970_2000" & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  select(c("cat.AI")) %>%
  table() %>% as.data.frame() %>%
  mutate(perc = round(Freq/sum(Freq)*100, 1),
         lab.y = rev(cumsum(rev(Freq))) - Freq*0.5) 

p <- ggplot(pie)+
    geom_bar(aes(x = "", y = Freq, fill = cat.AI), width = 1, stat = "identity")+
    coord_polar("y", start = 0)+
    scale_fill_manual(values = col.cat)+
    geom_label_repel(aes(x = "", y = lab.y, label = perc), nudge_x = 0.5)+
    labs(fill = "", title = "Multimodel CMIP6")+
    theme_void()+theme(legend.position = "right")

print(p)

sdmm <- cmip6 %>% subset(period == "1970_2000"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))


tab <- cmip6 %>% subset(period == "1970_2000"   & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)


knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	5.2	8.6	13.2	6.3	33.3	37.4	29.1	0.3
CESM	5.3	10.7	14.1	5.2	35.3	35.0	29.5	0.3
CMCC	5.1	7.4	10.1	5.4	28.0	42.2	29.5	0.3
CMCC-ESM2	6.4	7.9	7.7	2.5	24.5	24.4	50.8	0.3
CNRM	5.8	12.2	11.3	6.5	35.8	38.3	25.7	0.2
EC-Earth3	8.4	9.5	10.8	4.4	33.1	35.7	31.0	0.2
FGOALS	5.9	10.9	10.0	4.9	31.7	25.1	42.9	0.2
GFDL-ESM4	5.8	8.9	8.9	4.1	27.7	35.3	36.8	0.2
INM	4.4	8.4	13.3	5.7	31.8	40.1	27.8	0.3
INM-CM5	3.4	8.5	12.0	5.6	29.5	43.0	27.3	0.3
MPI	8.7	11.1	10.2	4.3	34.3	34.0	31.5	0.2
MRI	6.1	10.8	9.2	3.9	30.0	38.9	30.9	0.3
NorESM-2-MM	5.2	10.1	17.1	5.6	38.0	33.4	28.3	0.3
mean	5.8	9.6	11.4	5.0	31.8	35.6	32.4	0.3
sd	1.4	1.5	2.5	1.1	6.5	5.7	7.1	0.1

3.3 1985-2010

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "1985_2015" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "1985-2010", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "1985_2015"   & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))


tab <- cmip6 %>% subset(period == "1985_2015"   & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)


knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	5.3	8.4	13.8	6.6	34.1	37.2	28.3	0.3
CESM	4.8	10.8	14.8	5.4	35.8	35.4	28.5	0.3
CMCC	5.2	7.6	10.4	5.8	29.0	42.9	27.9	0.3
CMCC-ESM2	6.2	8.5	7.7	2.5	24.9	24.5	50.3	0.3
CNRM	5.8	12.1	11.3	6.2	35.4	39.2	25.2	0.2
EC-Earth3	8.3	9.8	10.7	4.4	33.2	36.7	29.9	0.2
FGOALS	5.8	11.4	10.4	4.9	32.5	25.3	42.0	0.2
GFDL-ESM4	5.9	9.1	9.1	3.9	28.0	36.5	35.3	0.2
INM	4.0	7.7	11.3	5.2	28.2	36.5	35.0	0.3
INM-CM5	3.5	8.3	12.0	5.4	29.2	40.6	29.9	0.3
MPI	8.7	11.1	10.3	4.4	34.5	34.1	31.0	0.2
MRI	6.2	10.7	9.1	3.8	29.8	38.4	31.6	0.3
NorESM-2-MM	4.7	9.9	15.9	5.6	36.1	30.3	33.2	0.3
mean	5.7	9.6	11.3	4.9	31.5	35.2	32.9	0.3
sd	1.5	1.5	2.3	1.1	6.4	5.5	6.7	0.1

3.4 SSP 2-4.5

3.4.1 2030-2060

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2030_2060" & model == "SSP245" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2030-2060", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2030_2060" & model == "SSP245" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2030_2060" & model == "SSP245"   & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	5.0	9.2	14.6	6.5	35.3	38.1	26.2	0.3
CESM	4.7	12.4	15.3	5.6	38.0	35.6	26.2	0.3
CMCC	4.4	8.4	11.0	6.1	29.9	44.3	25.5	0.3
CMCC-ESM2	6.1	9.1	8.0	2.8	26.0	25.3	48.4	0.3
CNRM	6.0	12.7	12.1	6.8	37.6	39.8	22.3	0.2
EC-Earth3	7.5	9.2	11.1	4.5	32.3	36.3	31.2	0.2
FGOALS	6.9	13.3	18.5	7.8	46.5	30.3	22.9	0.2
GFDL-ESM4	5.4	10.1	9.5	4.4	29.4	36.9	33.4	0.2
INM	4.6	8.4	13.8	5.9	32.7	40.4	26.7	0.3
INM-CM5	3.5	8.6	12.2	5.4	29.7	43.1	26.9	0.3
MPI	9.0	11.9	11.1	4.5	36.5	35.1	28.2	0.2
MRI	6.2	11.2	9.3	4.1	30.8	40.0	28.9	0.3
NorESM-2-MM	5.3	11.8	16.9	5.6	39.6	33.5	26.7	0.3
mean	4.0	7.3	8.7	3.7	23.7	25.6	48.0	2.6
sd	1.0	1.2	2.2	0.9	5.3	3.6	4.6	0.0

3.4.2 2070-2100

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2070_2100" & model == "SSP245" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2070-2100", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2070_2100" & model == "SSP245" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2070_2100" & model == "SSP245"   & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	4.8	10.5	15.0	7.3	37.6	38.2	24.1	0.3
CESM	6.1	11.6	15.4	5.4	38.5	36.9	24.3	0.3
CMCC	4.3	8.0	12.3	6.4	31.0	46.6	22.1	0.3
CMCC-ESM2	6.2	9.5	8.5	2.8	27.0	26.8	46.0	0.3
CNRM	6.3	12.9	12.4	6.8	38.4	41.0	20.4	0.2
EC-Earth3	7.6	9.1	11.0	4.4	32.1	37.3	30.3	0.2
FGOALS	6.9	14.0	18.1	8.1	47.1	30.3	22.3	0.2
GFDL-ESM4	5.9	9.7	9.8	4.4	29.8	38.4	31.6	0.2
INM	4.8	8.2	14.1	5.8	32.9	40.6	26.2	0.3
INM-CM5	3.7	9.5	11.2	5.9	30.3	43.5	25.9	0.3
MPI	9.0	11.9	11.1	4.5	36.5	35.1	28.2	0.2
MRI	6.5	11.0	9.6	4.0	31.1	40.8	27.9	0.3
NorESM-2-MM	5.9	11.7	17.0	5.9	40.5	34.0	25.3	0.3
mean	4.2	7.4	8.8	3.8	24.2	26.2	47.0	2.6
sd	1.0	1.3	2.1	1.0	5.4	3.7	4.5	0.0

3.5 SSP 3-7.0

3.5.1 2030-2060

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2030_2060" & model == "SSP370" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2030-2060", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2030_2060" & model == "SSP370" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2030_2060" & model == "SSP370"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	4.9	9.6	14.3	6.8	35.6	37.4	26.8	0.3
CESM	5.2	11.4	15.6	5.4	37.6	35.7	26.5	0.3
CMCC	4.2	8.3	11.3	5.9	29.7	44.1	25.9	0.3
CMCC-ESM2	3.3	8.4	9.0	3.6	24.3	26.9	48.6	0.3
CNRM	6.2	12.9	12.9	7.1	39.1	39.6	21.0	0.2
EC-Earth3	7.9	9.7	11.6	4.8	34.0	38.5	27.3	0.2
FGOALS	6.1	11.9	11.9	5.3	35.2	27.5	37.0	0.2
GFDL-ESM4	6.3	9.4	9.9	4.1	29.7	35.9	34.1	0.2
INM	4.6	8.3	13.5	5.9	32.3	40.3	27.1	0.3
INM-CM5	3.4	8.8	12.1	5.6	29.9	43.2	26.5	0.3
MPI	8.7	12.2	11.0	4.5	36.4	36.4	26.9	0.2
MRI	6.5	10.8	9.4	3.6	30.3	40.2	29.2	0.3
NorESM-2-MM	5.4	11.0	17.9	5.7	40.0	33.2	26.5	0.3
mean	3.9	7.1	8.6	3.7	23.3	25.6	48.5	2.6
sd	1.1	1.1	1.8	0.8	4.8	3.6	4.8	0.0

3.5.2 2070-2100

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2070_2100" & model == "SSP370" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2070-2100", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2070_2100" & model == "SSP370" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2070_2100" & model == "SSP370"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	5.5	10.5	15.4	6.9	38.3	37.3	24.2	0.3
CESM	6.0	12.2	15.6	5.7	39.5	37.4	22.9	0.3
CMCC	5.2	9.1	12.5	7.3	34.1	44.1	21.6	0.3
CMCC-ESM2	6.9	9.9	8.1	3.0	27.9	25.9	45.9	0.3
CNRM	7.3	12.9	13.2	6.7	40.1	39.3	20.3	0.2
EC-Earth3	7.6	9.8	12.4	4.7	34.5	40.3	24.9	0.2
FGOALS	6.7	13.3	13.3	4.4	37.7	28.0	34.1	0.2
GFDL-ESM4	7.5	8.5	10.9	4.7	31.6	37.0	31.1	0.2
INM	6.0	8.2	14.7	6.1	35.0	39.2	25.5	0.3
INM-CM5	3.9	10.3	14.1	6.1	34.4	40.4	24.9	0.3
MPI	8.7	12.2	11.0	4.5	36.4	36.4	26.9	0.2
MRI	7.2	11.9	9.9	4.1	33.1	39.0	27.6	0.3
NorESM-2-MM	5.9	11.8	17.5	5.9	41.1	33.7	24.9	0.3
mean	4.5	7.5	9.0	3.8	24.8	25.6	47.0	2.6
sd	0.9	1.2	1.8	0.9	4.8	3.5	4.7	0.0

3.6 SSP 5-8.5

3.6.1 2030-2060

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2030_2060" & model == "SSP585" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2030-2060", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2030_2060" & model == "SSP585" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2030_2060" & model == "SSP585"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(14,15), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	4.8	9.2	15.0	7.5	36.5	37.6	25.6	0.3
CESM	5.8	11.5	15.6	5.2	38.1	36.0	25.7	0.3
CMCC	4.4	8.0	11.5	6.5	30.4	44.8	24.5	0.3
CMCC-ESM2	6.5	9.3	8.0	2.5	26.3	25.4	48.0	0.3
CNRM	6.2	13.3	12.4	6.9	38.8	40.0	21.0	0.2
EC-Earth3	7.8	9.4	11.8	4.6	33.6	38.8	27.4	0.2
FGOALS	6.8	13.5	18.1	7.7	46.1	30.8	22.8	0.2
GFDL-ESM4	6.1	9.4	9.7	4.7	29.9	36.4	33.4	0.2
INM	4.6	8.4	14.0	6.0	33.0	40.6	26.2	0.3
INM-CM5	3.5	8.6	12.4	5.8	30.3	43.6	25.8	0.3
MPI	8.6	12.4	11.1	4.5	36.6	37.0	26.2	0.2
MRI	6.0	11.1	9.6	4.0	30.7	40.8	28.1	0.3
NorESM-2-MM	5.5	11.4	17.6	5.9	40.4	33.2	26.2	0.3
mean	4.1	7.3	8.9	3.8	24.1	26.0	47.3	2.6
sd	1.0	1.3	2.1	1.0	5.4	3.7	4.7	0.0

3.6.2 2070-2100

map_list <- list()

for(i in unique(cmip6$source)){
  g <- ggplot(data = subset(cmip6, period == "2070_2100" & model == "SSP585" & source == i)) + geom_raster(aes(x=lon, y = lat,  fill = cat.AI))+
    scale_fill_manual(values = col.cat, na.translate = F)+
    labs(title = paste(i, "2070-2100", sep = ", "), fill = "")+
    theme_void()
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 3, nrow = 5, common.legend = T, legend = "bottom")

sdmm <- cmip6 %>% subset(period == "2070_2100" & model == "SSP585" &!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  ungroup() %>% group_by(cat.AI) %>% summarise(mean = round(mean(percent, na.rm = T),1), sd = round(sd(percent, na.rm = T),1)) %>%
  t() %>% as.data.frame() %>% setNames(c(.[1,1:6], "NA")) %>% .[-1, ] %>% sapply(as.numeric) %>% as.data.frame() %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")]), "source" = c("mean","sd")) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

tab <- cmip6 %>% subset(period == "2070_2100" & model == "SSP585"  & !Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55) %>%
  group_by(source, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(source) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(source~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("source","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA")) %>% 
  rbind(sdmm)

knitr::kable(tab, digits = 2, escape = F)  %>%
  kable_styling(bootstrap_options = "bordered") %>% 
  column_spec(c(6,9), italic = T, include_thead = T) %>% row_spec(c(13,14), bold = T) 

source	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
CAS-ESM2	6.1	10.8	15.9	6.5	39.3	37.5	22.9	0.3
CMCC	4.4	9.3	13.5	6.7	33.9	46.7	19.2	0.3
CMCC-ESM2	5.6	10.5	8.9	2.9	27.9	28.1	43.8	0.3
CNRM	6.0	13.6	12.8	6.4	38.8	40.9	20.1	0.2
EC-Earth3	7.7	9.6	13.0	4.9	35.2	40.5	23.9	0.2
FGOALS	15.6	15.0	18.3	6.8	55.7	24.4	19.8	0.2
GFDL-ESM4	6.5	10.0	10.8	4.8	32.1	37.6	30.1	0.2
INM	4.4	8.7	14.8	6.0	33.9	41.7	24.1	0.3
INM-CM5	3.9	9.8	12.0	6.0	31.7	44.1	23.9	0.3
MPI	8.7	12.3	11.1	4.7	36.8	36.8	26.2	0.2
MRI	7.1	11.5	10.6	4.1	33.3	40.1	26.4	0.3
NorESM-2-MM	6.4	12.8	17.4	6.5	43.1	34.3	22.3	0.3
mean	4.8	7.8	9.2	3.8	25.6	26.2	45.5	2.6
sd	2.1	1.3	2.0	0.8	6.2	4.4	4.6	0.0

4 Agreement between models

4.1 Summary dataframe

4.1.1 Difference with AI ref for all SSP

test <- cmip6ext %>% subset(model %in% c("historical","SSP585") & period %in% c("1970_2000", "2070_2100") & lon == 25.5 & lat == 10.5) %>% select(c("source","lon","lat","model","period","AI","cat.AI","AI.ref"))
                                                                               
                                                                               
test2 <- cmip6s %>% subset(model == "SSP585" & period == "2070_2100" & lon == 25.5 & lat == 10.5) %>% 
  select(c("model","period","AI.mean","cat.AI","AI.ref.mean","diff.AI.mean", "nb.models", "diff.AI.mean"))

test2

cmip6 <- read.table("cmip6_AI.txt")

cmip6.ref <- cmip6 %>% filter(period == "1970_2000") %>% select(c("lon","lat","model","period","source","AI","spr","t")) %>% setNames(c("lon","lat","model","period","source","AI.ref","spr.ref","t.ref"))

cmip6ext <- merge(cmip6, select(cmip6.ref, c("lon","lat","source","AI.ref","spr.ref","t.ref")) , by = c("lon","lat","source"), all = T) %>%  mutate(diff.AI = AI - AI.ref, diff.pr = spr - spr.ref, diff.t = t - t.ref)

write.table(cmip6ext, "cmip6ext.txt")

cmip6 <- read.table("cmip6ext.txt")
cmip6s <- cmip6 %>% 
  group_by(lon,lat,Continent, Type, Name, Acronym, lm, model, period, z) %>% 
  dplyr::summarise(t.mean = mean(t, na.rm = T), t.sd = sd(t, na.rm = T), # t in °C
                   pr.mean = mean(spr, na.rm = T) , pr.sd = sd(spr, na.rm = T),  # pr in mm/y
                   sfcWind.mean = mean(sfcWind, na.rm = T), sfcWind.sd = sd(sfcWind, na.rm = T),
                   Rn.mean = mean(Rn, na.rm = T), Rn.sd = sd(Rn, na.rm = T), # in MJ/m2/day
                   ea.mean = mean(ea, na.rm = T), ea.sd = sd(ea, na.rm = T), # in kPa
                   ET0 = mean(sET0, na.rm = T), ET0.sd = sd(sET0, na.rm = T), # in mm/y
                   AI.mean = mean(AI, na.rm = T), AI.sd = sd(AI, na.rm = T),
                   AI.ref.mean = mean(AI.ref, na.rm = T), AI.ref.sd = sd(AI.ref, na.rm = T),
                   diff.AI.mean = mean(diff.AI, na.rm = T), diff.AI.sd = sd(diff.AI, na.rm = T),
                   spr.ref.mean = mean(spr.ref, na.rm = T), spr.ref.sd = sd(spr.ref, na.rm = T),
                   diff.pr.mean = mean(diff.pr, na.rm = T), diff.pr.sd = sd(diff.pr, na.rm = T),
                   t.ref.mean = mean(t.ref, na.rm = T), t.ref.sd = sd(t.ref, na.rm =T),
                   diff.t.mean = mean(diff.t, na.rm = T), diff.t.sd = sd(diff.t, na.rm = T),
                   Hyperarid = sum(cat.AI == "Hyperarid"),
                   Arid = sum(cat.AI == "Arid"),
                   'Semi-arid' = sum(cat.AI == "Semi-arid"),
                   'Dry subhumid' = sum(cat.AI == "Dry subhumid") ,
                   'Humid' = sum(cat.AI == "Humid"),
                   'Cold' = sum(cat.AI == "Cold"),
                   'NA' = sum(is.na(cat.AI)),
                   nb.models = n(),
                    plus.AI = table(diff.AI > 0)["TRUE"],
                    minus.AI = table(diff.AI < 0)["TRUE"]) %>%
ungroup() %>% as.data.frame()

cmip6s$cat.maj <- colnames(select(cmip6s, c("Hyperarid","Arid",'Semi-arid',"Dry subhumid",'Humid',"Cold","NA")))[max.col(select(cmip6s, c("Hyperarid","Arid",'Semi-arid',"Dry subhumid",'Humid',"Cold","NA")),ties.method = "first")]
cmip6s$nb.maj <- cmip6s %>% select(c("Hyperarid","Arid",'Semi-arid',"Dry subhumid",'Humid',"Cold","NA")) %>% apply(1, max)
cmip6s$prop.maj <- with(cmip6s, nb.maj / nb.models)



breaks.unesco <- c(-Inf,0,0.03,0.2,0.5,0.65,Inf)
cat.unesco <- c("(0,0.03]"= "Hyperarid", "(0.03,0.2]" = "Arid", "(0.2,0.5]" = "Semi-arid", "(0.5,0.65]" = "Dry subhumid","(0.65, Inf]" = "Humid", "(-Inf,0]" = "Cold")

cmip6s$cat.AI <- cmip6s$AI.mean %>% cut(breaks = breaks.unesco) %>% revalue(cat.unesco)
cmip6s$cat.AI[which(cmip6s$ET0 < 400)] <- "Cold"

cmip6s$cat.AI.ref <- cmip6s$AI.ref.mean %>% cut(breaks = breaks.unesco) %>% revalue(cat.unesco)
cmip6s$cat.AI.ref[which(cmip6s$ET0 < 400)] <- "Cold"

cmip6s$same.catAI <- with(cmip6s, ifelse(cat.maj == cat.AI, "same", "diff"))

cmip6s$cat.dryer <- with(cmip6s, 
          ifelse(cat.AI.ref == "Arid" & cat.AI %in% c("Cold","Humid","Dry subhumid","Semi-arid"), "wetter", 
          ifelse(cat.AI.ref == "Arid" & cat.AI == "Hyperarid", "dryer",
          ifelse(cat.AI.ref == "Semi-arid" & cat.AI %in% c("Cold","Humid","Dry subhumid"), "wetter",
          ifelse(cat.AI.ref == "Semi-arid" & cat.AI %in% c("Hyperarid","Arid"), "dryer",
          ifelse(cat.AI.ref == "Dry subhumid" & cat.AI %in% c("Cold","Humid"), "wetter",
          ifelse(cat.AI.ref == "Dry subhumid" & cat.AI %in% c("Semi-arid","Arid","Hyperarid"), "dryer",
          ifelse(cat.AI.ref == "Humid" & cat.AI == "Cold","wetter",
          ifelse(cat.AI.ref == "Humid" & cat.AI %in% c("Dry subhumid","Semi-arid","Arid","Hyperarid"), "dryer",
          ifelse(cat.AI.ref == "Cold" & cat.AI %in% c("Humid","Dry subhumid","Semi-arid","Arid","Hyperarid"), "dryer",
          ifelse(cat.AI.ref == "Hyperarid" & cat.AI %in% c("Cold","Humid","Dry subhumid","Semi-arid","Arid"), "wetter",
          ifelse(cat.AI.ref == cat.AI, "same",
          NA))))))))))))

cmip6sref <- cmip6s %>% filter(period == "1970_2000") %>% select(c("lon","lat","model","period","cat.maj")) %>% setNames(c("lon","lat","model","period","cat.maj.ref"))

cmip6s <- merge(cmip6s, select(cmip6sref, c("lon","lat", "model","period","cat.maj.ref")), 
                 by = c("lon","lat","model","period"), all = T) 

cmip6s$cat.maj.dryer <- with(cmip6s, 
        ifelse(cat.maj.ref == "Arid" & cat.maj %in% c("Cold","Humid","Dry subhumid","Semi-arid"), "wetter", 
        ifelse(cat.maj.ref == "Arid" & cat.maj == "Hyperarid", "dryer",
        ifelse(cat.maj.ref == "Semi-arid" & cat.maj %in% c("Cold","Humid","Dry subhumid"), "wetter",
        ifelse(cat.maj.ref == "Semi-arid" & cat.maj %in% c("Hyperarid","Arid"), "dryer",
        ifelse(cat.maj.ref == "Dry subhumid" & cat.maj %in% c("Cold","Humid"), "wetter",
        ifelse(cat.maj.ref == "Dry subhumid" & cat.maj %in% c("Semi-arid","Arid","Hyperarid"), "dryer",
        ifelse(cat.maj.ref == "Humid" & cat.maj == "Cold","wetter",
        ifelse(cat.maj.ref == "Humid" & cat.maj %in% c("Dry subhumid","Semi-arid","Arid","Hyperarid"), "dryer",
        ifelse(cat.maj.ref == "Cold" & cat.maj %in% c("Humid","Dry subhumid","Semi-arid","Arid","Hyperarid"), "dryer",
        ifelse(cat.maj.ref == "Hyperarid" & cat.maj %in% c("Cold","Humid","Dry subhumid","Semi-arid","Arid"), "wetter",
        ifelse(cat.maj.ref == cat.maj, "same",
        NA))))))))))))


write.table(cmip6s, "cmip6s.txt")

cmip6s <- read.table("cmip6s.txt")

4.2 Compare cat mean and cat maj

cat mean = category of aridity obtained by using the multimodel mean AI cat maj = category of aridity obtained by using the category majoritarily obtained in the 13 models

4.2.1 Where are cat mean and cat maj different

cmip6s %>% group_by(model, period) %>% summarise(maj.equal.mean = table(same.catAI)[2], maj.diff.mean = table(same.catAI)[1], prop.same = round(maj.equal.mean/n()*100,1)) %>% knitr::kable(col.names = c("Model","Period","Number of gridcells in which cat.maj = cat.mean", "Number of gridcells in which cat.maj /= cat.mean", "Proportion of cells with the same category (%)")) %>% kableExtra::kable_styling(bootstrap_options = "bordered")

Model	Period	Number of gridcells in which cat.maj = cat.mean	Number of gridcells in which cat.maj /= cat.mean	Proportion of cells with the same category (%)
SSP245	2030_2060	19884	1798	89.3
SSP245	2070_2100	19890	1792	89.3
SSP370	2030_2060	19988	1694	89.8
SSP370	2070_2100	19911	1771	89.4
SSP585	2030_2060	19799	1883	88.9
SSP585	2070_2100	19796	1886	88.9
historical	1850_1880	19907	1775	89.4
historical	1970_2000	19880	1802	89.3
historical	1985_2015	19936	1746	89.6

cmip6s$cat.AI <- factor(cmip6s$cat.AI, levels =  c("Hyperarid","Arid","Semi-arid","Dry subhumid", "Humid","Cold"))
cmip6s$cat.maj <- factor(cmip6s$cat.maj, levels =  c("Hyperarid","Arid","Semi-arid","Dry subhumid", "Humid","Cold"))

ggplot(subset(cmip6s, cat.AI != "Cold"))+
  geom_col(aes(x= period, y = cat.AI, fill = cat.AI), just = -0.2, width = 0.3)+
  geom_col(aes(x = period, y = cat.maj, fill = cat.maj), just = 1.2, width = 0.3)+
  scale_fill_manual(values = col.cat, na.translate = F)+
  labs(fill = "", y ="")+
  facet_grid(rows = vars(model), switch = "y")+
  theme_minimal()+theme(axis.text.y = element_blank())

tab <- rbind(table(cmip6s$cat.AI)/length(cmip6s$cat.AI)*100, table(cmip6s$cat.maj)/length(cmip6s$cat.maj)*100) %>% t() %>% as.data.frame() %>% setNames(c("Cat mean", "Cat maj"))
knitr::kable(tab, digits =2, caption = "Proportion of gridcells in each aridity category") %>% kable_styling(bootstrap_options = "bordered")

Proportion of gridcells in each aridity category

	Cat mean	Cat maj
Hyperarid	3.51	4.18
Arid	6.25	7.43
Semi-arid	7.63	9.19
Dry subhumid	3.39	0.84
Humid	28.75	28.53
Cold	47.90	47.23

map_list <- list()

for(i in unique(cmip6s$model)){
  df <- subset(cmip6s, model == i)
  for(j in unique(df$period)){
    index <- paste(i, j, sep = "")
    g <- ggplot()+geom_raster(data = subset(df, period == j & same.catAI == "diff"), aes(x = lon, y = lat, fill = cat.AI))+
      borders(colour = "grey60")+
      scale_fill_manual(values =  col.cat, na.translate = F)+
      labs(title = paste(i, j , sep = " "), fill = "")+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }
}

col1 <- ggarrange(plotlist = map_list, nrow = 9, ncol = 1, common.legend = T, legend = "bottom") %>% annotate_figure(top = "cat.mean for gridcells in which\ncat.maj and cat.mean are different")

map_list2 <- list()

for(i in unique(cmip6s$model)){
  df <- subset(cmip6s, model == i)
  for(j in unique(df$period)){
    index <- paste(i, j, sep = "")
    g <- ggplot()+geom_raster(data = subset(df, period == j & same.catAI == "diff"), aes(x = lon, y = lat, fill = cat.maj))+
      borders(colour = "grey60")+
      scale_fill_manual(values =  col.cat, na.translate = F)+
      labs(title = paste(i, j , sep = " "), fill = "")+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list2[[index]] <- g
  }
}

col2 <- ggarrange(plotlist = map_list2, nrow = 9, ncol = 1, common.legend = T, legend = "bottom") %>% annotate_figure(top = "cat.maj for gridcells in which\ncat.maj and cat.mean are different")

ggarrange(col1,col2,ncol = 2, common.legend = T)

4.2.2 Plot cat maj and cat mean for reference period

ggpubr::ggarrange(plotlist = list(
  
  ggplot() + 
    geom_raster(data = subset(cmip6s, period == "1970_2000"), aes(x=lon, y = lat,  fill = cat.maj))+
    geom_point(data = subset(cmip6s, period == "1970_2000" & nb.maj > 10), aes(x =lon, y = lat), shape = "'", size = 0.5, col = "grey60")+
    scale_fill_manual(values = col.cat)+
    labs(title = "Category chosen in majority by models", fill = "")+
    ylim(-55,90)+
    theme_void(base_size = 15)+
    theme(legend.position = "none"),
  
  ggplot() + 
    geom_raster(data = subset(cmip6s, period == "1970_2000"), aes(x=lon, y = lat,  fill = cat.AI))+
   geom_point(data = subset(cmip6s, period == "1970_2000" & nb.maj > 10), aes(x =lon, y = lat), shape = "'", size = 0.5, col = "grey60")+
    scale_fill_manual(values = col.cat)+
    labs(title = "Multimodel average categories of AI", fill = "")+
    ylim(-55,90)+
    theme_void(base_size = 15)+
    theme(legend.position = "none")),
  
  ncol = 2)

map_list <- list()

for(i in c("1850_1880","1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = cat.maj))+
    geom_point(data = subset(cmip6s, period == i & model == "historical" & nb.maj > 10), aes(x= lon, y = lat), shape = "'", size = 0.5, col = "grey60")+
    scale_fill_manual(values =  col.cat, na.translate = F)+
    labs(title = i, fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")
  map_list[[i]] <- g
}


map_list2 <- list()

for(i in c("1850_1880","1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = cat.AI))+
   # geom_point(data = subset(cmip6s, period == i & model == "historical" & nb.maj > 10), aes(x= lon, y = lat), size = 0.5, col = "grey60")+
    scale_fill_manual(values =  col.cat, na.translate = F)+
    labs(title = i, fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")
  map_list2[[i]] <- g
}


ggarrange(
  annotate_figure(ggarrange(plotlist = map_list, nrow = 3, common.legend = T, legend = "bottom"), top = "Majority category among models"),
  annotate_figure(ggarrange(plotlist = map_list2, nrow = 3, common.legend = T, legend = "bottom"), top = "Multimodel average category"),
  ncol = 2, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = cat.maj))+
         geom_point(data = subset(cmip6s, period == j & model == i & nb.maj > 10), aes(x= lon, y = lat), shape = "'", size = 0.5, col = "grey60")+
      scale_fill_manual(values = col.cat, na.translate = F)+
      labs(fill = "Cat AI maj", title = paste(i, j, sep = ", "))+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}


map_list2 <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = cat.AI))+
    #     geom_point(data = subset(cmip6s, period == j & model == i & nb.maj > 10), aes(x= lon, y = lat), shape = "'", size = 0.5, col = "grey60")+
      scale_fill_manual(values = col.cat, na.translate = F)+
      labs(fill = "Cat AI mean", title = paste(i, j, sep = ", "))+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list2[[index]] <- g
  }}

ggarrange(plotlist = map_list, nrow = 3, ncol = 2, common.legend = T, legend = "bottom")

ggarrange(plotlist = map_list2, nrow = 3, ncol = 2, common.legend = T, legend = "bottom")

4.3 Multimodel maps of cat mean

cmip6s$cat.AI <- factor(cmip6s$cat.AI, levels = c("Hyperarid","Arid","Semi-arid","Dry subhumid","Humid","Cold"))

map_list <- list()

for(i in c("1850_1880", "1970_2000","1985_2015")){
  g <- ggplot() + 
    geom_tile(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = cat.AI))+
    geom_point(data = subset(cmip6s,nb.maj > 10), aes(x = lon, y = lat), shape = ".", col = "grey20", size = 0.1)+
    borders(colour = "grey60")+
    scale_fill_manual(values =  col.cat, na.translate = F)+
    labs(title = paste("Category AI",i , sep = " "), fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")
  map_list[[i]] <- g
}

ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + 
      geom_tile(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = cat.AI))+
      geom_point(data = subset(cmip6s,nb.maj > 10), aes(x = lon, y = lat), shape = ".", col = "grey20", size = 0.1)+
      borders(colour = "grey60")+
      scale_fill_manual(values =  col.cat, na.translate = F)+
      labs(title = paste("Category AI",index, sep = " "), fill = "")+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

5 Maps multimodel AI

quantile(cmip6s$AI.mean, probs = seq(0,1,0.1), na.rm = T) 
##           0%          10%          20%          30%          40%          50% 
## 2.902197e-03 1.993698e-01 5.922059e-01 1.019621e+00 1.457766e+00 2.127218e+00 
##          60%          70%          80%          90%         100% 
## 3.578980e+00 1.168470e+01 2.438867e+01 4.856852e+01 5.329875e+06

ai.breaks <- c(0,0.03,0.2,0.5,0.65,1,10,20,30)
colscale <- colorvec[-c(5,6)]

map_list <- list()

for(i in c("1850_1880", "1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = AI.mean))+
    borders(colour = "grey60")+
    binned_scale(aesthetics = "fill", breaks = ai.breaks, palette = function(x) colscale,
                 guide = guide_legend(label.theme = element_text(angle = 0)))+
    labs(title = i, fill = "Aridity index")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right")
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = AI.mean))+
      borders(colour = "grey60")+
      binned_scale(aesthetics = "fill", breaks = ai.breaks, palette = function(x) colscale,
                   guide = guide_legend(label.theme = element_text(angle = 0)))+
      labs(fill = "Aridity index", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

6 AI changes over time and scenarios

v<-NULL
v2<-NULL
for (N in 5:30) {
  v<-c(v,sum(sapply(floor(2*N/3):N,function(i) choose(N,i)*(1/2)^N*2)))# proba d'avoir 2/3 de même signe à partir de N modèles, lorsque N change
  v2<-c(v2,sum(sapply(floor(2*N/3):N,function(i) choose(N,i)*(1/2)^N)))# proba d'avoir 2/3 de signe >0 à partir de N modèles, lorsque N change
}


plot(5:30,v,type="l",ylab="Proba(>2/3 models...)",xlab="Nb of models")
lines(5:30,v2,col="red")

6.1 Changes AI for summary CMIP6

6.2 Boxplot of changes AI

ggpubr::ggarrange(plotlist = list(
  ggplot(subset(cmip6s, diff.AI.mean < 1 & diff.AI.mean > -1 & model == "historical"))+
    geom_boxplot(aes(x=period, y = - diff.AI.mean, col = period), outliers = F)+
    scale_color_manual(values = c("#1A9850","#66BD63","#A6D96A"))+
    labs(x = "", y = "Difference of AI between each period and\nthe reference period 1970-2000")+
    facet_grid(cols = vars(model))+
    theme_minimal() + theme(legend.position = "none", axis.text.x = element_text(angle = 45, vjust = 1, hjust=1)),
  ggplot(subset(cmip6s, diff.AI.mean < 1 & diff.AI.mean > -1 & model != "historical"))+
    geom_boxplot(aes(x=period, y = diff.AI.mean, col = period), outliers = F)+
    scale_color_manual(values = c("#FDAE61","#F46D43"))+
    scale_y_continuous(position = "right")+
    labs(y = "", x = "")+
    facet_grid(cols = vars(model))+
    theme_minimal() + theme(legend.position = "none", axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))),
  widths = c(1.5,3),
  ncol = 2, nrow = 1
)

6.3 Strictly wetter or dryer (reference period: 1970-2000)

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
  index <- paste(i, j, sep = ", ")
  g <- ggplot()+
  geom_point(data = subset(cmip6s, model == i & period == j & plus.AI > 10), aes(x = lon, y = lat, col = "plus"), shape = "+")+
  geom_point(data = subset(cmip6s, model == i & period == j & minus.AI > 10), aes(x = lon, y = lat, col = "minus"), shape = "-")+
    scale_color_manual(values = colorvec[c(1,11)])+
    labs(color = "", x = "", y = "", title = index)+
  theme_minimal()
  map_list[[index]] <- g
  }
}

ggarrange(plotlist = map_list, ncol = 2, nrow = 3, common.legend = T)

map_list <- list()

for(i in c("1850_1880","1985_2015")){
    index <- paste(i, j, sep = " , ")
  g <- ggplot() + geom_tile(data = subset(cmip6s, period == i & model == "historical" & diff.AI.mean > 0), aes(x=lon, y = lat,  fill = "wetter"))+
    geom_tile(data = subset(cmip6s, period == i & model == "historical" & diff.AI.mean < 0), aes(x=lon, y = lat,  fill = "dryer"))+
    geom_point(data = subset(cmip6s, period == i & model == "historical" & plus.AI > 10), aes(x=lon, y = lat), shape = "+", size = 0.5)+
     geom_point(data = subset(cmip6s, period == i & model == "historical" & minus.AI > 10), aes(x=lon, y = lat), shape = "-", size = 0.5)+
    borders(colour = "grey60")+
    scale_fill_manual(values =  c("#F46D43", "#74ADD1"), na.translate = F)+
    labs(title = paste("AI changes between",i , "and 1970-2000", sep = " "), fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")
  map_list[[i]] <- g
  }

ggarrange(plotlist = map_list, ncol = 2)

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i & diff.AI.mean > 0), aes(x=lon, y = lat,  fill = "wetter"))+
      geom_tile(data = subset(cmip6s, period == j & model == i & diff.AI.mean < 0), aes(x=lon, y = lat,  fill = "dryer"))+
       geom_point(data = subset(cmip6s, period == j & model == i & plus.AI > 10), aes(x=lon, y = lat), shape = "+", col = "grey20", size = 0.5)+
       geom_point(data = subset(cmip6s, period == j & model == i & minus.AI > 10), aes(x=lon, y = lat), shape = "-", col = "grey20", size = 0.5)+
      borders(colour = "grey60")+
      scale_fill_manual(values = c("#F46D43", "#74ADD1"), na.translate = F)+
      labs(fill = "Compared to 1970-2000", title = i)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}



ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

6.4 Wetter or dryer, binned (reference period: 1970-2000)

ano.quantiles <- cmip6s %>% group_by(model, period) %>% summarise(q10 = quantile(diff.AI.mean, probs = 0.1, na.rm =T), 
                                                                  q25 = quantile(diff.AI.mean, probs = 0.25, na.rm =T), 
                                                                  q50 = quantile(diff.AI.mean, probs = 0.5, na.rm =T), 
                                                                  q75 = quantile(diff.AI.mean, probs = 0.75, na.rm =T), 
                                                                  q90 = quantile(diff.AI.mean, probs = 0.9, na.rm =T)) 

ano.quantiles <- cmip6s %>% summarise(q10 = quantile(diff.AI.mean, probs = 0.1, na.rm =T), 
                                      q25 = quantile(diff.AI.mean, probs = 0.25, na.rm =T), 
                                      q50 = quantile(diff.AI.mean, probs = 0.5, na.rm =T), 
                                      q75 = quantile(diff.AI.mean, probs = 0.75, na.rm =T), 
                                      q90 = quantile(diff.AI.mean, probs = 0.9, na.rm =T)) 
qbreaks <- c(-5, -1,-0.1, -0.01, 0, 0.01, 0.1, 1,5)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#fddbc7", "#f4a582", "#d6604d", "#b2182b", "#67001f")

“diff.AI.mean” is the multimodel AI average minus the reference AI (period 1970-2000). More blue: AI mean is superior to AI ref, hence wetter than the ref. More red: the anomalie is negative, hence AI mean is inferior to AI ref, hence dryer.

AS OF NOW THE SAME QUANTILES ARE USED FOR ALL MAPS.

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  #quant = subset(ano.quantiles, model == "historical" & period == i)
  #b = quant[1,c(3:7)] %>% as.numeric
  b = qbreaks
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = diff.AI.mean))+
    borders(colour = "grey60")+
           geom_point(data = subset(cmip6s, period == i & model == "historical" & plus.AI > 10), aes(x=lon, y = lat), shape = "+", col = "grey20", size = 0.5)+
       geom_point(data = subset(cmip6s,  period == i & model == "historical"  & minus.AI > 10), aes(x=lon, y = lat), shape = "-", col = "grey20", size = 0.5)+
    binned_scale(aesthetics = "fill", breaks = b, palette = function(x) rev(colscale),
                 guide = guide_legend(label.theme = element_text(angle = 0)))+
    labs(title = i, fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right")
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 2, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    b = qbreaks
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = diff.AI.mean))+
      borders(colour = "grey60")+
                 geom_point(data = subset(cmip6s, period == i & model == "historical" & plus.AI > 10), aes(x=lon, y = lat), shape = "+", col = "grey20", size = 0.5)+
       geom_point(data = subset(cmip6s,  period == i & model == "historical"  & minus.AI > 10), aes(x=lon, y = lat), shape = "-", col = "grey20", size = 0.5)+
      binned_scale(aesthetics = "fill", breaks = b, palette = function(x) rev(colscale),
                   guide = guide_legend(label.theme = element_text(angle = 0)))+
      labs(fill = "Compared to 1970-2000", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}



ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

6.5 Wetter or dryer, binned, % (reference period: 1970-2000)

cmip6s$diff.AI.mean.percent <- with(cmip6s, (AI.mean-AI.ref.mean)/AI.ref.mean*100)


pbreaks <- c(-40, -30,-20, -10, 0, 10, 20, 30, 40)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef","white","white", "#fddbc7", "#f4a582", "#d6604d", "#b2182b")

“diff.AI.mean” is the multimodel AI average minus the reference AI (period 1970-2000). More blue: AI mean is superior to AI ref, hence wetter than the ref. More red: the anomalie is negative, hence AI mean is inferior to AI ref, hence dryer.

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  #quant = subset(ano.quantiles, model == "historical" & period == i)
  #b = quant[1,c(3:7)] %>% as.numeric
  b = pbreaks
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = diff.AI.mean.percent))+
    borders(colour = "grey60")+
           geom_point(data = subset(cmip6s, period == i & model == "historical" & plus.AI > 10), aes(x=lon, y = lat), shape = "+", col = "grey20", size = 0.5)+
       geom_point(data = subset(cmip6s,  period == i & model == "historical"  & minus.AI > 10), aes(x=lon, y = lat), shape = "-", col = "grey20", size = 0.5)+
    binned_scale(aesthetics = "fill", breaks = b, palette = function(x) rev(colscale))+
    labs(title = i, fill = "")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right")
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 2, common.legend = T, legend = "bottom") %>% annotate_figure(top = "% change of AI index compared to 1970-2000")

#,guide = guide_legend(label.theme = element_text(angle = 0))

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    b = pbreaks
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = diff.AI.mean.percent))+
      borders(colour = "grey60", size = 0.5)+
                 geom_point(data = subset(cmip6s, period == j & model == i & plus.AI > 10), aes(x=lon, y = lat), shape = "+", col = "grey20", size = 0.5)+
       geom_point(data = subset(cmip6s,  period == j & model == i & minus.AI > 10), aes(x=lon, y = lat), shape = "-", col = "grey20", size = 0.5)+
      binned_scale(aesthetics = "fill", breaks = b, palette = function(x) rev(colscale))+
      labs(fill = "% change of AI index\ncompared to 1970-2000", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}


ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

7 Difference between future and reference category of AI, using cat maj

7.1 Detailed changes and colors cat maj

cmip6s$change.catmaj <- with(cmip6s, ifelse(cat.maj == cat.maj.ref, 0, 1))

cmip6s$change.cat.maj <- with(cmip6s, paste(cat.maj.ref, cat.maj, sep = " to "))

cmip6s$change.cat.maj <- with(cmip6s, ifelse(cat.maj.ref == "Arid" & cat.maj %in% c("Cold","Humid"),"Arid to cold or humid",
         ifelse(cat.maj.ref == "Arid" & cat.maj %in% c("Semi-Arid","Dry subhumid"), "Arid to semi-arid or dry subhumid",
         ifelse(cat.maj.ref == "Arid" & cat.maj == "Hyperarid", "Arid to hyperarid",
         ifelse(cat.maj.ref == "Semi-arid" & cat.maj == "Dry subhumid", "Semi-arid to dry subhumid",
         ifelse(cat.maj.ref == "Semi-arid" & cat.maj == "Arid", "Semi-arid to arid",
         ifelse(cat.maj.ref == "Semi-arid" & cat.maj == "Hyperarid", "Semi-arid to hyperarid",
         ifelse(cat.maj.ref == "Semi-arid" & cat.maj %in% c("Humid","Cold"), "Semi-arid to humid or cold",
         ifelse(cat.maj.ref == "Dry subhumid" & cat.maj %in% c("Semi-arid","Arid","Hyperarid"), "Dry subhumid to arid, semi-arid or hyperarid",
         ifelse(cat.maj.ref == "Dry subhumid" & cat.maj %in% c("Humid","Cold"), "Dry subhumid to humid or cold",
         ifelse(cat.maj.ref == "Humid" & cat.maj == "Cold","Humid to cold",
         ifelse(cat.maj.ref == "Humid" & cat.maj %in% c("Dry subhumid", "Semi-arid"), "Humid to dry-subhumid or semi-arid",
         ifelse(cat.maj.ref == "Humid" & cat.maj %in% c("Arid", "Hyperarid"), "Humid to arid or hyperarid",
         ifelse(cat.maj.ref == "Cold" & cat.maj == "Humid", "Cold to humid", 
         ifelse(cat.maj.ref == "Cold" & cat.maj %in% c("Dry subhumid","Semi-arid","Arid","Hyperarid"), "Cold to dryland",
         ifelse(cat.maj.ref == "Hyperarid" & cat.maj == "Arid", "Hyperarid to arid",
         ifelse(cat.maj.ref == "Hyperarid" & cat.maj %in% c("Semi-arid","Dry-subhumid"),"Hyperarid to semi arid or dry subhumid",
         ifelse(cat.maj.ref == "Hyperarid" & cat.maj %in% c("Humid", "Cold"), "Hyperarid to humid or cold", NA))))))))))))))))))

cmip6s$change.cat.maj <- factor(cmip6s$change.cat.maj, 
                                levels = c("Arid to hyperarid","Semi-arid to arid", "Dry subhumid to arid, semi-arid or hyperarid", "Humid to dry-subhumid or semi-arid", "Cold to humid","Hyperarid to arid", "Arid to semi-arid or dry subhumid", "Dry subhumid to humid or cold", "Semi-arid to dry subhumid", "Semi-arid to humid or cold"))
table(cmip6s$change.cat.maj)
## 
##                            Arid to hyperarid 
##                                            0 
##                            Semi-arid to arid 
##                                            0 
## Dry subhumid to arid, semi-arid or hyperarid 
##                                            0 
##           Humid to dry-subhumid or semi-arid 
##                                            0 
##                                Cold to humid 
##                                            0 
##                            Hyperarid to arid 
##                                            0 
##            Arid to semi-arid or dry subhumid 
##                                            0 
##                Dry subhumid to humid or cold 
##                                            0 
##                    Semi-arid to dry subhumid 
##                                            0 
##                   Semi-arid to humid or cold 
##                                            0

colors_changes <- c("Hyperarid to arid" = "#E7D7B6", 
                    "Arid to semi-arid or dry subhumid" = colorvec[7],
                    "Dry subhumid to humid or cold" = colorvec[8],
                    "Semi-arid to dry subhumid" = colorvec[9],
                    "Semi-arid to humid or cold" = colorvec[10],
                    "Humid to cold" = colorvec[11],
                    "Cold to humid" = "#E7D7B6", 
                    "Humid to dry-subhumid or semi-arid" = colorvec[5],
                    "Dry subhumid to arid, semi-arid or hyperarid" = colorvec[4],
                    "Semi-arid to arid" = colorvec[3],
                    "Arid to hyperarid" = colorvec[1])

7.2 Changes towards dryer categories cat maj

7.2.1 Historical

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  g <- ggplot() + 
    geom_raster(data = subset(cmip6s, period == i & model == "historical" & change.catmaj == 1 & cat.maj.dryer == "dryer"), 
                aes(x=lon, y = lat,  fill = change.cat.maj))+
    borders(colour = "grey60")+
    scale_fill_manual(values =  colors_changes, na.translate = F)+
    labs(title = i, fill = "Change of category\ncompared to\n1970-2000")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
  map_list[[i]] <- g
}


ggarrange(plotlist = map_list, ncol = 2, common.legend = T, legend = "bottom")

7.2.1.1 Future

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i & change.catmaj == 1 & cat.maj.dryer == "dryer"), 
                              aes(x=lon, y = lat,  fill = change.cat.maj))+
      #   geom_raster(data = subset(cmip6s, period == i & model == "historical" & change.catmaj == 1 & nb.maj > 10), aes(x= lon, y = lat, fill = "'"))+
      borders(colour = "grey60")+
      scale_fill_manual(values = colors_changes, na.translate = F)+
      #  scale_fill_viridis_d()+
      labs(fill = "Change of category\ncompared to\n1970-2000", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
    map_list[[index]] <- g
  }}

ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

7.3 Changes towards wetter categories cat maj

7.3.0.1 Historical

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical" & change.catmaj == 1 & cat.maj.dryer == "wetter"), 
                              aes(x=lon, y = lat,  fill = change.cat.maj))+
    borders(colour = "grey60", ylim = c(-40,40))+
    scale_fill_manual(values =  colors_changes, na.translate = F)+
    labs(title = i, fill = "Change of category\ncompared to\n1970-2000")+
    theme_void()+
    ylim(-50,50)+
    theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
  map_list[[i]] <- g
}

ggarrange(plotlist = map_list, ncol = 2, common.legend = T, legend = "bottom")

7.3.0.2 Future

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i & change.catmaj == 1 & cat.maj.dryer == "wetter"), 
                              aes(x=lon, y = lat,  fill = change.cat.maj))+
      borders(colour = "grey60", ylim = c(-40,40))+
      scale_fill_manual(values = colors_changes, na.translate = F)+
      labs(fill = "Change of category\ncompared to\n1970-2000", title = index)+
      theme_void()+ylim(-50,50)+
      theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list, nrow = 3, ncol = 2, common.legend = T, legend = "bottom")

8 Difference between future and reference category of AI, using cat mean

8.1 Detailed changes and colors

cmip6s$change.catmean <- with(cmip6s, ifelse(cat.AI == cat.AI.ref, 0, 1))

cmip6s$change.cat.mean <- with(cmip6s, paste(cat.AI.ref, cat.AI, sep = " to "))
unique(cmip6s$change.cat)
## NULL

# cmip6s$change.cat.mean <- with(cmip6s, 
#       ifelse(cat.AI.ref == "Arid" & cat.AI %in% c("Cold","Humid"),"Arid to cold or humid",
#       ifelse(cat.AI.ref == "Arid" & cat.AI %in% c("Semi-Arid","Dry subhumid"), "Arid to semi-arid or dry subhumid",
#       ifelse(cat.AI.ref == "Arid" & cat.AI == "Hyperarid", "Arid to hyperarid",
#       ifelse(cat.AI.ref == "Semi-arid" & cat.AI == "Dry subhumid", "Semi-arid to dry subhumid",
#       ifelse(cat.AI.ref == "Semi-arid" & cat.AI == "Arid", "Semi-arid to arid",
#       ifelse(cat.AI.ref == "Semi-arid" & cat.AI == "Hyperarid", "Semi-arid to hyperarid",
#       ifelse(cat.AI.ref == "Semi-arid" & cat.AI %in% c("Humid","Cold"), "Semi-arid to humid or cold",
#       ifelse(cat.AI.ref == "Dry subhumid" & cat.AI %in% c("Semi-arid","Arid","Hyperarid"), "Dry subhumid to arid, semi-arid or hyperarid",
#       ifelse(cat.AI.ref == "Dry subhumid" & cat.AI %in% c("Humid","Cold"), "Dry subhumid to humid or cold",
#       ifelse(cat.AI.ref == "Humid" & cat.AI == "Cold","Humid to cold",
#       ifelse(cat.AI.ref == "Humid" & cat.AI %in% c("Dry subhumid","Semi-arid"), "Humid to dry-subhumid or semi-arid",
#       ifelse(cat.AI.ref == "Humid" & cat.AI %in% c("Arid", "Hyperarid"), "Humid to arid or hyperarid",
#       ifelse(cat.AI.ref == "Cold" & cat.AI == "Humid", "Cold to humid", 
#       ifelse(cat.AI.ref == "Cold" & cat.AI %in% c("Dry subhumid","Semi-arid","Arid","Hyperarid"), "Cold to dryland",
#       ifelse(cat.AI.ref == "Hyperarid" & cat.AI == "Arid", "Hyperarid to arid",
#       ifelse(cat.AI.ref == "Hyperarid" & cat.AI %in% c("Semi-arid","Dry-subhumid"),"Hyperarid to semi arid or dry subhumid",
#       ifelse(cat.AI.ref == "Hyperarid" & cat.AI %in% c("Humid", "Cold"), "Hyperarid to humid or cold", 
#       NA))))))))))))))))))

table(cmip6s$change.cat.mean)
## 
##                 Arid to Arid            Arid to Hyperarid 
##                        11080                          639 
##            Arid to Semi-arid                 Cold to Cold 
##                          340                        95966 
## Dry subhumid to Dry subhumid        Dry subhumid to Humid 
##                         4022                          288 
##    Dry subhumid to Semi-arid        Humid to Dry subhumid 
##                         2242                         2484 
##               Humid to Humid           Humid to Semi-arid 
##                        57318                          198 
##            Hyperarid to Arid       Hyperarid to Hyperarid 
##                          240                         6402 
##                     NA to NA            Semi-arid to Arid 
##                         5130                         1209 
##    Semi-arid to Dry subhumid       Semi-arid to Semi-arid 
##                          277                        12514
cmip6s$change.cat.mean <- factor(cmip6s$change.cat.mean, 
                                  levels = c("Arid to Hyperarid","Semi-arid to Arid", "Dry subhumid to Semi-arid", "Humid to Semi-arid", "Humid to Dry subhumid", "Cold to Humid","Hyperarid to arid", "Arid to Semi-arid", "Dry subhumid to Humid", "Semi-arid to Dry subhumid", "Semi-arid to Humid"))

colors_changes <- c("Hyperarid to Arid" = "#E7D7B6", 
                     "Arid to Semi-arid" = colorvec[7],
                     "Dry subhumid to Humid" = colorvec[8],
                     "Semi-arid to Dry subhumid" = colorvec[9],
                     "Semi-arid to Humid" = colorvec[10],
                     "Humid to Cold" = colorvec[11],
                     "Cold to Humid" = "#E7D7B6",
                     "Humid to Dry subhumid" = colorvec[5],
                     "Humid to Semi-arid" = colorvec[4],
                     "Dry subhumid to Semi-arid" = colorvec[3],
                     "Semi-arid to Arid" = colorvec[2],
                     "Arid to Hyperarid" = colorvec[1])

8.2 Changes towards dryer categories cat mean

8.2.1 Historical

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  g <- ggplot() + 
    geom_raster(data = subset(cmip6s, period == i & model == "historical" & change.catmean == 1 & cat.dryer == "dryer"), 
                aes(x=lon, y = lat,  fill = change.cat.mean))+
    borders(colour = "grey60")+
    scale_fill_manual(values =  colors_changes, na.translate = F)+
    labs(title = i, fill = "Change of category\ncompared to\n1970-2000")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
  map_list[[i]] <- g
}

ggarrange(plotlist = map_list, ncol = 2, common.legend = T, legend = "bottom")

8.2.1.1 Future

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i & change.catmean == 1 & cat.dryer == "dryer"), 
                              aes(x=lon, y = lat,  fill = change.cat.mean))+
      borders(colour = "grey60")+
      scale_fill_manual(values = colors_changes)+
      labs(fill = "Change of category\ncompared to\n1970-2000", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
    map_list[[index]] <- g
  }}

ggarrange(plotlist = map_list, ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

8.3 Changes towards wetter categories cat mean

8.3.1 Historical

map_list <- list()

for(i in c("1850_1880","1985_2015")){
  g <- ggplot() + geom_tile(data = subset(cmip6s, period == i & model == "historical" & change.catmaj == 1 & cat.dryer == "wetter"), 
                            aes(x=lon, y = lat,  fill = change.cat.mean))+
    borders(colour = "grey60", ylim = c(-40,40))+
    scale_fill_manual(values =  colors_changes, na.translate = F)+
    labs(title = i, fill = "Change of category\ncompared to\n1970-2000")+
    theme_void()+ylim(-50,50)+
    theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
  map_list[[i]] <- g
}

ggarrange(plotlist = map_list, ncol = 2, common.legend = T, legend = "bottom")

8.3.2 Future

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i & change.catmean == 1 & cat.dryer == "wetter"), 
                              aes(x=lon, y = lat,  fill = change.cat.mean))+
      borders(colour = "grey60", ylim = c(-40,40))+
      scale_fill_manual(values = colors_changes, na.translate = F)+
      labs(fill = "Change of category\ncompared to\n1970-2000", title = index)+
      theme_void()+
      theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 2))
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

8.4 Big map of SSP585 2070-2100

g <- ggplot() + 
  geom_tile(data = subset(cmip6s, period == "2070_2100" & model == "SSP585" & change.catmean == 1), aes(x=lon, y = lat,  fill = change.cat.mean))+
      borders(colour = "grey60")+
      scale_fill_manual(values = colors_changes, na.translate = F)+
      labs(fill = "Change of aridity category\nfor the period 2070-2100\nand SSP 5-8.5,\ncompared to 1970-2000", title = "")+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")+guides(fill = guide_legend(nrow = 3))
   
print(g)

9 Tables

9.1 Changes of aridity category using cat mean

tab <- cmip6s %>% group_by(model, period, Continent) %>% summarise(Cold = table(cat.AI)[1], Hyperarid = table(cat.AI)[2], Arid = table(cat.AI)[3], Semiarid = table(cat.AI)[4], Drysubhumid = table(cat.AI)[5], Humid = table(cat.AI)[6], total = n()) 
write.table(tab, "tab.catAI.txt")

tab.future <- cmip6s %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(period, model, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(period, model) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(period + model ~cat.AI) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

kable(tab.future[order(tab.future$model, decreasing = F),]) %>% kable_styling(bootstrap_options = "bordered") %>%
  column_spec(c(7,10), italic = T, include_thead = T) %>% row_spec(2, bold = T) 

	model	period	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
1	historical	1850_1880	3.4	6.0	7.1	3.3	19.8	28.8	48.9	2.6
2	historical	1970_2000	3.3	6.0	7.0	3.3	19.6	27.8	49.9	2.6
3	historical	1985_2015	3.3	6.1	6.9	3.3	19.6	27.3	50.6	2.6
4	SSP245	2030_2060	3.5	6.3	7.8	3.3	20.9	28.8	47.7	2.6
7	SSP245	2070_2100	3.6	6.4	7.8	3.5	21.3	29.4	46.6	2.6
5	SSP370	2030_2060	3.2	6.7	7.6	3.0	20.5	28.8	48.1	2.6
8	SSP370	2070_2100	4.0	6.3	8.3	3.7	22.3	28.5	46.6	2.6
6	SSP585	2030_2060	3.6	6.2	7.9	3.5	21.2	29.3	47.0	2.6
9	SSP585	2070_2100	3.8	6.4	8.4	3.6	22.2	30.1	45.1	2.6

# ggplot(tab.future, aes(x = period))+
#   geom_point(aes(y = Hyperarid, col = "Hyperarid", shape = model))+
#   geom_line(aes(y=Hyperarid, col = "Hyperarid", group = model, lty = model))+
#   geom_point(aes(y = Arid, col = "Arid", shape = model))+
#   geom_line(aes(y=Arid, col = "Arid", group = model, lty = model))+
#   geom_point(aes(y = get('Semi-arid'), col = "Semi-arid", shape = model))+
#   geom_line(aes(y=get('Semi-arid'), col = "Semi-arid", group = model, lty = model))+
#   geom_point(aes(y = get('Dry subhumid'), col = "Dry subhumid", shape = model))+
#   geom_line(aes(y=get('Dry subhumid'), col = "Dry subhumid", group = model, lty = model))+
#   geom_point(aes(y = Humid, col = "Humid", shape = model))+
#   geom_line(aes(y= Humid, col = "Humid", group = model, lty = model))+
#   scale_color_manual(values = col.cat)+
#   theme_minimal()

cmip6$cat.AI <- factor(cmip6$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))
tab.percent <- cmip6 %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, cat.AI, source) %>%
  summarise(count = n()) %>%
  ungroup() %>% 
  group_by(period, model, cat.AI) %>% 
  summarise(mmmean = mean(count, na.rm = T), mmsd = sd(count, na.rm = T)) %>%
  mutate(percent = round(mmmean/sum(mmmean)*100, 1), sd.percent = round(mmsd/sum(mmmean)*100,1)) 
write.table(tab.percent, "tab.percent.txt")

tab.percent <- read.table("tab.percent.txt")

df245 <- rbind(subset(tab.percent, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.percent, model == "SSP245")) %>%
  group_by(period, model) %>% 
  mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

df245$cat.AI <- factor(df245$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))

b245 <- ggplot(data = df245, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = paste(percent, " ± ", sd.percent, sep = "")), size = 3)+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 2-4.5", y = "%", x = "", col = "", fill = "")+
  theme_minimal()

df370 <- rbind(subset(tab.percent, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.percent, model == "SSP370")) %>%
  group_by(period, model) %>% 
  mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

df370$cat.AI <- factor(df370$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))

b370 <- ggplot(data = df370, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = paste(percent, " ± ", sd.percent, sep = "")), size = 3)+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 3-7.0", y = "%", x = "", col = "", fill = "")+
  theme_minimal()


df585 <- rbind(subset(tab.percent, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.percent, model == "SSP585")) %>%
  group_by(period, model) %>% 
  mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

df585$cat.AI <- factor(df585$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))

b585 <- ggplot(data = df585, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = paste(percent, " ± ", sd.percent, sep = "")), size = 3)+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 5-8.5", y = "%", x = "", col = "", fill = "")+
  theme_minimal()


ggarrange(plotlist = list(b245, b370, b585), common.legend = T, legend = "bottom", ncol = 3)

9.2 Changes of aridity category using cat maj

tab.maj <- cmip6s %>% group_by(model, period, Continent) %>% summarise(Cold = table(cat.maj)[1], Hyperarid = table(cat.maj)[2], Arid = table(cat.maj)[3], Semiarid = table(cat.maj)[4], Drysubhumid = table(cat.maj)[5], Humid = table(cat.maj)[6], total = n()) 
write.table(tab.maj, "tab.catmaj.txt")

tab.future.maj <- cmip6s %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC")) %>%
  group_by(period, model, cat.maj) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(period, model) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL) %>%
  reshape2::dcast(period + model ~cat.maj) %>%
  mutate("Sum drylands" = rowSums(.[c("Hyperarid","Arid","Semi-arid","Dry subhumid")])) %>%
  select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))

kable(tab.future.maj[order(tab.future.maj$model, decreasing = F),]) %>% kable_styling(bootstrap_options = "bordered") %>%
  column_spec(c(7,10), italic = T, include_thead = T) %>% row_spec(2, bold = T) 

	model	period	Hyperarid	Arid	Semi-arid	Dry subhumid	Sum drylands	Humid	Cold	NA
1	historical	1850_1880	4.2	7.1	8.7	0.8	20.8	28.3	48.3	2.6
2	historical	1970_2000	4.1	7.1	8.6	0.7	20.5	28.1	48.9	2.6
3	historical	1985_2015	4.0	7.0	8.4	0.9	20.3	27.3	49.9	2.6
4	SSP245	2030_2060	4.0	7.8	9.3	0.8	21.9	28.4	47.2	2.6
7	SSP245	2070_2100	4.2	7.7	9.4	0.9	22.2	29.1	46.1	2.6
5	SSP370	2030_2060	4.1	7.2	9.4	0.8	21.5	28.6	47.3	2.6
8	SSP370	2070_2100	4.5	7.7	10.0	0.9	23.1	28.4	45.9	2.6
6	SSP585	2030_2060	4.2	7.5	9.6	1.0	22.3	28.8	46.3	2.6
9	SSP585	2070_2100	4.5	7.9	9.4	0.9	22.7	30.0	44.6	2.6

tab.flow.catmaj <- cmip6s %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, cat.maj) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(period, model) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL)

df245 <- rbind(subset(tab.flow.catmaj, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow.catmaj, model == "SSP245")) %>%
  mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.maj %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b245 <- ggplot(data = df245, aes(x = period, y = percent))+
  geom_col(aes(group = cat.maj, col = cat.maj, fill = cat.maj))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 2-4.5", y = "%", x = "")+
  theme_minimal()

df370 <- rbind(subset(tab.flow.catmaj, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow.catmaj, model == "SSP370")) %>% mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.maj %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b370 <- ggplot(data = df370, aes(x = period, y = percent))+
  geom_col(aes(group = cat.maj, col = cat.maj, fill = cat.maj))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 3-7.0", y = "%", x = "")+
  theme_minimal()


df585 <- rbind(subset(tab.flow.catmaj, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow.catmaj, model == "SSP585")) %>% mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.maj %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b585 <- ggplot(data = df585, aes(x = period, y = percent))+
  geom_col(aes(group = cat.maj, col = cat.maj, fill = cat.maj))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  labs(title = "SSP 5-8.5", y = "%", x = "")+
  theme_minimal()


ggarrange(plotlist = list(b245, b370, b585), common.legend = T, legend = "bottom", ncol = 3)

9.3 Changes by continent, cat mean

cmip6$cat.AI <- factor(cmip6$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))

tab.percent.cont <- cmip6 %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, Continent, cat.AI, source) %>%
  summarise(count = n()) %>%
  ungroup() %>% 
  group_by(period, model, Continent, cat.AI) %>% 
  summarise(mmmean = mean(count, na.rm = T), mmsd = sd(count, na.rm = T)) %>%
  mutate(percent = round(mmmean/sum(mmmean)*100, 1), sd.percent = round(mmsd/sum(mmmean)*100,1)) 
write.table(tab.percent.cont, "tab.percent.continent.txt")

tab.percent.cont <- read.table("tab.percent.continent.txt")

tab_list_percent <- list()
tab_list_sd <- list()

for(i in unique(tab.percent.cont$Continent)){
  tab.i <- subset(tab.percent.cont, Continent == i)
  
  df.percent <- tab.i %>% reshape2::dcast(period + model ~ cat.AI, value.var = "percent")
  
  df.sd <- tab.i %>% reshape2::dcast(period + model ~ cat.AI, value.var = "sd.percent")
  # sd for a sum: square root of sum of squared sd
  
  drycats <- which(names(df.percent) %in% c("Hyperarid","Arid","Semi-arid","Dry subhumid"))
  
  df.percent <- df.percent %>% mutate("Sum drylands" = rowSums(.[drycats], na.rm = T))
  df.sd <- df.sd %>% mutate("Sum drylands" = round(sqrt(rowSums((.[drycats])^2, na.rm = T)),2))
  
  missing <- setdiff(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"), names(df.percent))
  
  df.percent[, missing] <- 0
  df.sd[,missing] <- 0
  
  df.percent <- df.percent %>% select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))
  df.percent[is.na(df.percent)] <- 0
  
  tab_list_percent[[i]] <- df.percent
  
  df.sd <- df.sd %>% select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))
  df.sd[is.na(df.sd)] <- 0
  
  tab_list_sd[[i]] <- df.sd
}

k_list <- list()


kdiff_list <- list()

for(i in names(tab_list_percent)){
  
  df <- tab_list_percent[[i]] %>% select(c("model", "period"))
  
  
  
  prop.ref <- tab_list_percent[[i]] %>% filter(period == "1970_2000")
  prop.ref <- prop.ref %>% rbind(prop.ref[rep(1,dim(df)[1]-1),])
  df2 <- tab_list_percent[[i]][,c(3:10)] - prop.ref[,c(3:10)] %>% round(1)
  diff.prop.percent <- cbind(df, df2)
  
  sd.ref <- tab_list_sd[[i]] %>% filter(period == "1970_2000")
  sd.ref <- sd.ref %>% rbind(sd.ref[rep(1,dim(df)[1]-1),])
  df3 <- sqrt(tab_list_sd[[i]][,c(3:10)]^2 + sd.ref[,c(3:10)]^2) %>% round(1)
  diff.sd.percent <- cbind(df, df3)
  
  df$Hyperarid <- paste(tab_list_percent[[i]]$Hyperarid, tab_list_sd[[i]]$Hyperarid, sep = " ± ")
  df$Arid <- paste(tab_list_percent[[i]]$Arid, tab_list_sd[[i]]$Arid, sep = " ± ")
  df$'Semi-Arid' <- paste(tab_list_percent[[i]][,5], tab_list_sd[[i]][,5], sep = " ± ")
  df$'Dry subhumid' <- paste(tab_list_percent[[i]][,6], tab_list_sd[[i]][,6], sep = " ± ")
  df$'Sum drylands' <- paste(tab_list_percent[[i]][,7], tab_list_sd[[i]][,7], sep = " ± ")
  df$Humid <- paste(tab_list_percent[[i]]$Humid, tab_list_sd[[i]]$Humid, sep = " ± ")
  df$Cold <- paste(tab_list_percent[[i]]$Cold, tab_list_sd[[i]]$Cold, sep = " ± ")
  
  df.diff <- df
  df.diff$Hyperarid <- paste(round(diff.prop.percent$Hyperarid, 1), round(diff.sd.percent$Hyperarid, 1), sep = " ± ")
  df.diff$Arid <- paste(round(diff.prop.percent$Arid, 1), round(diff.sd.percent$Arid, 1), sep = " ± ")
  df.diff$'Semi-Arid' <- paste(round(diff.prop.percent[,5], 1), round(diff.sd.percent[,5], 1), sep = " ± ")
  df.diff$'Dry subhumid' <- paste(round(diff.prop.percent[,6], 1), round(diff.sd.percent[,6], 1), sep = " ± ")
  df.diff$'Sum drylands' <- paste(round(diff.prop.percent[,7],1), round(diff.sd.percent[,7], 1), sep = " ± ")
  df.diff$Humid <- paste(round(diff.prop.percent$Humid, 1), round(diff.sd.percent$Humid, 1), sep = " ± ")
  df.diff$Cold <- paste(round(diff.prop.percent$Cold, 1), round(diff.sd.percent$Cold, 1), sep = " ± ")
  
  k_list[[i]] <- kable(df[order(df$model, decreasing = F),], caption = i) %>% kable_styling(bootstrap_options = "bordered") %>%
    column_spec(8, italic = T, background = colorvec[5]) %>% row_spec(2, bold = T) 
  kdiff_list[[i]] <- kable(df.diff[order(df.diff$model, decreasing = F),], caption = i) %>% kable_styling(bootstrap_options = "bordered") %>% 
    column_spec(8, italic = T, background = colorvec[5]) %>% row_spec(2, bold = T) 
  
}

9.3.1 Proportion of aridity categories by continent

k_list

$AFRICA

AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	27.2 ± 4.7	17.1 ± 4.9	15.6 ± 3.6	6.3 ± 1.4	66.2 ± 7.81	33.7 ± 7.1	0 ± 0
2	historical	1970_2000	26.9 ± 4.7	16.7 ± 4.2	15.1 ± 3.4	6.2 ± 1.5	64.9 ± 7.32	32.8 ± 6.2	2.3 ± 0
3	historical	1985_2015	26.5 ± 4.8	17.1 ± 4.6	15.4 ± 3.3	6 ± 1.5	65 ± 7.57	32.7 ± 6.3	2.2 ± 0
4	SSP245	2030_2060	27 ± 4.8	18 ± 4.9	16.5 ± 3.5	6.4 ± 1.7	67.9 ± 7.89	32 ± 7.6	0 ± 0
7	SSP245	2070_2100	27.9 ± 4.7	17.4 ± 5	17.2 ± 3.4	6.2 ± 1.5	68.7 ± 7.8	31.3 ± 7.9	0 ± 0
5	SSP370	2030_2060	26 ± 5	16.6 ± 3.6	17.2 ± 1.8	6.8 ± 1	66.6 ± 6.5	32.3 ± 6.7	1.1 ± 0
8	SSP370	2070_2100	28 ± 3.5	17.6 ± 4.8	16.6 ± 3.1	6.2 ± 1.6	68.4 ± 6.89	30.4 ± 6.4	1.2 ± 0
6	SSP585	2030_2060	27.1 ± 4.6	17.6 ± 4.8	16.8 ± 3.5	6.3 ± 1.5	67.8 ± 7.66	32.1 ± 7.5	0 ± 0
9	SSP585	2070_2100	29.7 ± 10.5	17.3 ± 5.3	16.2 ± 4.2	6 ± 1.5	69.2 ± 12.58	30.9 ± 8.6	0 ± 0

$ASIA

ASIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	3.9 ± 2.3	12.2 ± 1.9	9.2 ± 1.5	4 ± 1.4	29.3 ± 3.62	32.4 ± 7.7	38 ± 9.2
2	historical	1970_2000	3.8 ± 2.1	12 ± 1.8	8.9 ± 1.5	3.8 ± 1.1	28.5 ± 3.33	30.7 ± 8.5	40.5 ± 9.6
3	historical	1985_2015	3.7 ± 2.2	11.8 ± 1.6	8.6 ± 1.5	3.7 ± 1.1	27.8 ± 3.3	30.1 ± 7.6	41.8 ± 8.8
4	SSP245	2030_2060	3.5 ± 1.8	12.5 ± 1.6	9.4 ± 1.6	4.3 ± 1.7	29.7 ± 3.35	34.3 ± 8.1	35.8 ± 9.2
7	SSP245	2070_2100	3.7 ± 2	12.5 ± 1.7	9.5 ± 1.7	4.2 ± 1.7	29.9 ± 3.56	36.2 ± 8.6	33.7 ± 9.3
5	SSP370	2030_2060	3.4 ± 2.3	12.4 ± 1.5	8.9 ± 1.7	3.7 ± 1.3	28.4 ± 3.48	34 ± 7.2	37.3 ± 9.5
8	SSP370	2070_2100	4.6 ± 1.9	12.5 ± 1.5	9.6 ± 1.5	3.8 ± 1.2	30.5 ± 3.09	35.2 ± 9.2	34.1 ± 9.7
6	SSP585	2030_2060	3.8 ± 2	12.5 ± 1.6	9.5 ± 1.8	4.2 ± 1.6	30 ± 3.52	35.4 ± 8.4	34.4 ± 9.3
9	SSP585	2070_2100	4.8 ± 3.2	13.1 ± 3.2	10.4 ± 3.4	4.2 ± 1.5	32.5 ± 5.86	36.5 ± 11.7	30.7 ± 9.9

$CENTRAL-AMERICA

CENTRAL-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	8.4 ± 3.6	30.1 ± 5.1	16.1 ± 8.7	54.6 ± 10.71	45 ± 6.1	0 ± 0
2	historical	1970_2000	0 ± 0	8.4 ± 4.1	29.7 ± 7.6	14.7 ± 7.8	52.8 ± 11.64	44.9 ± 12	2 ± 0
3	historical	1985_2015	0 ± 0	9 ± 4.2	30.7 ± 7	14.8 ± 7.4	54.5 ± 11.02	43.8 ± 11.5	1.4 ± 0
4	SSP245	2030_2060	0 ± 0	11.9 ± 6.2	36.3 ± 8	13.9 ± 6.2	62.1 ± 11.87	37.6 ± 8.6	0 ± 0
7	SSP245	2070_2100	0 ± 0	12.3 ± 5.5	38.1 ± 8.5	14.4 ± 7.9	64.8 ± 12.84	34.8 ± 8.6	0 ± 0
5	SSP370	2030_2060	0 ± 0	12.6 ± 6.3	31.4 ± 7.1	15 ± 7.9	59 ± 12.35	40.7 ± 11.8	0 ± 0
8	SSP370	2070_2100	4 ± 5.2	16 ± 5.5	35.5 ± 9.5	13.2 ± 5.2	68.7 ± 13.21	30.9 ± 10.6	0 ± 0
6	SSP585	2030_2060	0 ± 0	13 ± 7	37.6 ± 9	14.1 ± 7.1	64.7 ± 13.43	35 ± 9.4	0 ± 0
9	SSP585	2070_2100	2 ± 2.4	16.6 ± 8	41.2 ± 7.7	13.3 ± 7.7	73.1 ± 13.72	26.6 ± 8.8	0 ± 0

$EUROPE

EUROPE

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	1.9 ± 1.2	9.6 ± 4.1	4.5 ± 2.8	16 ± 5.11	49 ± 15.3	34.8 ± 20.2
2	historical	1970_2000	0 ± 0	1.5 ± 1	8.9 ± 3.4	3.9 ± 2.1	14.3 ± 4.12	44.9 ± 16.1	40.6 ± 20.4
3	historical	1985_2015	0 ± 0	1.3 ± 1	8.4 ± 2.9	3.8 ± 2	13.5 ± 3.66	46 ± 16.1	40.3 ± 19.9
4	SSP245	2030_2060	0 ± 0	1.9 ± 0.9	10.3 ± 3.6	5.1 ± 1.9	17.3 ± 4.17	52.7 ± 15.4	29.8 ± 19.7
7	SSP245	2070_2100	0 ± 0	2.2 ± 1.1	10.1 ± 4.8	5.1 ± 2.3	17.4 ± 5.44	55.4 ± 14.4	27 ± 18.6
5	SSP370	2030_2060	0 ± 0	2.1 ± 1	9.7 ± 4.7	4.8 ± 2.5	16.6 ± 5.42	52.7 ± 15.7	30.4 ± 20
8	SSP370	2070_2100	0 ± 0	3.1 ± 1.1	11.4 ± 5.9	6.4 ± 3.8	20.9 ± 7.1	52.8 ± 14.3	26 ± 18.6
6	SSP585	2030_2060	0 ± 0	2.3 ± 1.1	10.5 ± 3.7	4.9 ± 2.4	17.7 ± 4.55	55.2 ± 15.8	26.9 ± 19.9
9	SSP585	2070_2100	0 ± 0	3 ± 1.7	11.1 ± 5.4	6.1 ± 2.9	20.2 ± 6.36	57.2 ± 12.8	22.5 ± 17.3

$EUROPE-AFRICA

EUROPE-AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	5.9 ± 4.9	23.6 ± 7	21.3 ± 8.8	8.8 ± 3.3	59.6 ± 12.7	24.9 ± 14.2	14.5 ± 0
2	historical	1970_2000	7 ± 5.5	26.5 ± 7.5	23.3 ± 10.3	9.9 ± 3.2	66.7 ± 14.24	27.5 ± 16.5	4.5 ± 7.5
3	historical	1985_2015	7.1 ± 5.9	26.3 ± 7.1	22.2 ± 7.7	10.3 ± 4.2	65.9 ± 12.73	27.8 ± 14.3	5.2 ± 7
4	SSP245	2030_2060	9.9 ± 6.6	28.9 ± 8.6	26.7 ± 8.9	9.4 ± 3.4	74.9 ± 14.43	23.1 ± 16.2	0.9 ± 0
7	SSP245	2070_2100	11.7 ± 5.7	28.6 ± 8.8	27.2 ± 8.1	9.5 ± 3.1	77 ± 13.61	21.9 ± 13.9	0 ± 0
5	SSP370	2030_2060	10.3 ± 6.3	28.2 ± 8.4	26 ± 8.7	9 ± 2.8	73.5 ± 13.92	24.8 ± 14	0.5 ± 0
8	SSP370	2070_2100	13.2 ± 6.2	27.3 ± 8.6	28.3 ± 9.7	8.7 ± 3.6	77.5 ± 14.81	21.2 ± 16.2	0 ± 0
6	SSP585	2030_2060	11.1 ± 6.4	28.8 ± 8.3	27.5 ± 8.8	8.7 ± 3	76.1 ± 14.01	22.6 ± 14.1	0.2 ± 0
9	SSP585	2070_2100	15.2 ± 5.7	29.1 ± 10.1	27.9 ± 8.8	8 ± 2.8	80.2 ± 14.82	18.6 ± 13.6	0 ± 0

$NORTH-AMERICA

NORTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	0.4 ± 0.2	5.8 ± 5.6	4.3 ± 2.6	10.5 ± 6.18	35.1 ± 9	54.2 ± 10.9
2	historical	1970_2000	0 ± 0	0.4 ± 0.2	4.3 ± 4.4	4.2 ± 2.4	8.9 ± 5.02	34.5 ± 9.8	56.5 ± 11.7
3	historical	1985_2015	0 ± 0	0.3 ± 0.2	4.6 ± 3.8	4.1 ± 2.2	9 ± 4.4	33.6 ± 10.3	57.2 ± 11.5
4	SSP245	2030_2060	0 ± 0	0.4 ± 0.3	6.5 ± 5.5	4.7 ± 2.6	11.6 ± 6.09	37.5 ± 8.7	50.8 ± 11.1
7	SSP245	2070_2100	0 ± 0	0.5 ± 0.3	6.4 ± 5.1	5.1 ± 2.9	12 ± 5.87	39.7 ± 9	48.1 ± 11
5	SSP370	2030_2060	0 ± 0	0.5 ± 0.4	5.9 ± 4.8	4.6 ± 2.7	11 ± 5.52	37 ± 9.7	51.8 ± 12.2
8	SSP370	2070_2100	0.1 ± 0	0.8 ± 0.7	7.6 ± 4	5.2 ± 2	13.7 ± 4.53	38.4 ± 10.6	47.7 ± 11.2
6	SSP585	2030_2060	0 ± 0	0.5 ± 0.4	6.6 ± 5.2	5.2 ± 2.8	12.3 ± 5.92	38.2 ± 9.1	49.3 ± 11.7
9	SSP585	2070_2100	0.2 ± 0	1 ± 1.4	7 ± 6.1	5.2 ± 2.5	13.4 ± 6.74	42.1 ± 9.3	44.5 ± 11.3

$OCEANIA

OCEANIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	5.3 ± 0	32 ± 26	36.7 ± 17.1	6.6 ± 3.1	80.6 ± 31.27	18.9 ± 12.5	0.3 ± 0.2
2	historical	1970_2000	5.2 ± 0	31.3 ± 26.4	36.8 ± 17.4	6.8 ± 3.1	80.1 ± 31.77	19.3 ± 12.8	0.5 ± 0.4
3	historical	1985_2015	5.2 ± 0	33.4 ± 26.1	35.5 ± 17.9	6.8 ± 3.2	80.9 ± 31.81	18.4 ± 12.3	0.6 ± 0.3
4	SSP245	2030_2060	10.4 ± 0	35.8 ± 23.5	31.7 ± 15.1	5.7 ± 2.6	83.6 ± 28.05	15.9 ± 11.9	0.3 ± 0.3
7	SSP245	2070_2100	10.5 ± 0	37.9 ± 22.2	29.1 ± 13.1	6.2 ± 3.3	83.7 ± 25.99	16.2 ± 12	0.1 ± 0
5	SSP370	2030_2060	10.7 ± 0	34.8 ± 22.6	33.4 ± 15.3	5.8 ± 2.5	84.7 ± 27.41	14.9 ± 8.1	0.3 ± 0.2
8	SSP370	2070_2100	10.5 ± 10.2	38.1 ± 19.7	28.5 ± 13.8	5.7 ± 3.2	82.8 ± 26.32	16.8 ± 16.3	0.3 ± 0.3
6	SSP585	2030_2060	5.6 ± 3.6	35.6 ± 25.3	33.8 ± 16.7	6.7 ± 3.6	81.7 ± 30.74	17.7 ± 13.8	0.3 ± 0.3
9	SSP585	2070_2100	9.3 ± 6.4	39.4 ± 21.7	29.4 ± 12.9	5.5 ± 2.5	83.6 ± 26.16	16.2 ± 12.5	0 ± 0

$SOUTH-AMERICA

SOUTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0.1 ± 0.1	2.1 ± 1.7	15.2 ± 8.7	7.4 ± 2.6	24.8 ± 9.24	73.2 ± 11.7	1.9 ± 1
2	historical	1970_2000	0.1 ± 0.1	1.8 ± 1.3	13.9 ± 4.3	7.3 ± 2.3	23.1 ± 5.05	74.5 ± 7.3	2.3 ± 2.6
3	historical	1985_2015	0.2 ± 0.1	2 ± 1.2	14.5 ± 5.1	7.5 ± 2.4	24.2 ± 5.76	73.6 ± 7.7	2.1 ± 2.6
4	SSP245	2030_2060	0.2 ± 0.2	2.7 ± 2.2	17.4 ± 9.6	8.8 ± 3.6	29.1 ± 10.49	69.4 ± 13.1	1.3 ± 0.9
7	SSP245	2070_2100	0.2 ± 0.1	3.3 ± 3.5	18.3 ± 9.1	9.3 ± 4.6	31.1 ± 10.78	67.5 ± 14.7	1.3 ± 0.7
5	SSP370	2030_2060	0.2 ± 0.1	2.2 ± 1.7	16.3 ± 5.5	8.8 ± 3.7	27.5 ± 6.84	70.5 ± 8.9	1.9 ± 1.6
8	SSP370	2070_2100	0.2 ± 0.1	3.3 ± 2.5	18.6 ± 7	8.7 ± 3.8	30.8 ± 8.35	67.3 ± 10.8	1.8 ± 1.8
6	SSP585	2030_2060	0.2 ± 0.1	3.1 ± 2.5	18.6 ± 9	9.2 ± 4.4	31.1 ± 10.33	67.6 ± 13.1	1.2 ± 0.8
9	SSP585	2070_2100	0.6 ± 1	4.4 ± 3.8	19.6 ± 9.5	9.7 ± 3.9	34.3 ± 11	64.4 ± 12.6	1.1 ± 0.7

9.3.2 Changes in proportion of aridity categories by continent

kdiff_list

$AFRICA

AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0.3 ± 6.6	0.4 ± 6.5	0.5 ± 5	0.1 ± 2.1	1.3 ± 10.7	0.9 ± 9.4	-2.3 ± 0
2	historical	1970_2000	0 ± 6.6	0 ± 5.9	0 ± 4.8	0 ± 2.1	0 ± 10.4	0 ± 8.8	0 ± 0
3	historical	1985_2015	-0.4 ± 6.7	0.4 ± 6.2	0.3 ± 4.7	-0.2 ± 2.1	0.1 ± 10.5	-0.1 ± 8.8	-0.1 ± 0
4	SSP245	2030_2060	0.1 ± 6.7	1.3 ± 6.5	1.4 ± 4.9	0.2 ± 2.3	3 ± 10.8	-0.8 ± 9.8	-2.3 ± 0
7	SSP245	2070_2100	1 ± 6.6	0.7 ± 6.5	2.1 ± 4.8	0 ± 2.1	3.8 ± 10.7	-1.5 ± 10	-2.3 ± 0
5	SSP370	2030_2060	-0.9 ± 6.9	-0.1 ± 5.5	2.1 ± 3.8	0.6 ± 1.8	1.7 ± 9.8	-0.5 ± 9.1	-1.2 ± 0
8	SSP370	2070_2100	1.1 ± 5.9	0.9 ± 6.4	1.5 ± 4.6	0 ± 2.2	3.5 ± 10.1	-2.4 ± 8.9	-1.1 ± 0
6	SSP585	2030_2060	0.2 ± 6.6	0.9 ± 6.4	1.7 ± 4.9	0.1 ± 2.1	2.9 ± 10.6	-0.7 ± 9.7	-2.3 ± 0
9	SSP585	2070_2100	2.8 ± 11.5	0.6 ± 6.8	1.1 ± 5.4	-0.2 ± 2.1	4.3 ± 14.6	-1.9 ± 10.6	-2.3 ± 0

$ASIA

ASIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0.1 ± 3.1	0.2 ± 2.6	0.3 ± 2.1	0.2 ± 1.8	0.8 ± 4.9	1.7 ± 11.5	-2.5 ± 13.3
2	historical	1970_2000	0 ± 3	0 ± 2.5	0 ± 2.1	0 ± 1.6	0 ± 4.7	0 ± 12	0 ± 13.6
3	historical	1985_2015	-0.1 ± 3	-0.2 ± 2.4	-0.3 ± 2.1	-0.1 ± 1.6	-0.7 ± 4.7	-0.6 ± 11.4	1.3 ± 13
4	SSP245	2030_2060	-0.3 ± 2.8	0.5 ± 2.4	0.5 ± 2.2	0.5 ± 2	1.2 ± 4.7	3.6 ± 11.7	-4.7 ± 13.3
7	SSP245	2070_2100	-0.1 ± 2.9	0.5 ± 2.5	0.6 ± 2.3	0.4 ± 2	1.4 ± 4.9	5.5 ± 12.1	-6.8 ± 13.4
5	SSP370	2030_2060	-0.4 ± 3.1	0.4 ± 2.3	0 ± 2.3	-0.1 ± 1.7	-0.1 ± 4.8	3.3 ± 11.1	-3.2 ± 13.5
8	SSP370	2070_2100	0.8 ± 2.8	0.5 ± 2.3	0.7 ± 2.1	0 ± 1.6	2 ± 4.5	4.5 ± 12.5	-6.4 ± 13.6
6	SSP585	2030_2060	0 ± 2.9	0.5 ± 2.4	0.6 ± 2.3	0.4 ± 1.9	1.5 ± 4.8	4.7 ± 12	-6.1 ± 13.4
9	SSP585	2070_2100	1 ± 3.8	1.1 ± 3.7	1.5 ± 3.7	0.4 ± 1.9	4 ± 6.7	5.8 ± 14.5	-9.8 ± 13.8

$CENTRAL-AMERICA

CENTRAL-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	0 ± 5.5	0.4 ± 9.2	1.4 ± 11.7	1.8 ± 15.8	0.1 ± 13.5	-2 ± 0
2	historical	1970_2000	0 ± 0	0 ± 5.8	0 ± 10.7	0 ± 11	0 ± 16.5	0 ± 17	0 ± 0
3	historical	1985_2015	0 ± 0	0.6 ± 5.9	1 ± 10.3	0.1 ± 10.8	1.7 ± 16	-1.1 ± 16.6	-0.6 ± 0
4	SSP245	2030_2060	0 ± 0	3.5 ± 7.4	6.6 ± 11	-0.8 ± 10	9.3 ± 16.6	-7.3 ± 14.8	-2 ± 0
7	SSP245	2070_2100	0 ± 0	3.9 ± 6.9	8.4 ± 11.4	-0.3 ± 11.1	12 ± 17.3	-10.1 ± 14.8	-2 ± 0
5	SSP370	2030_2060	0 ± 0	4.2 ± 7.5	1.7 ± 10.4	0.3 ± 11.1	6.2 ± 17	-4.2 ± 16.8	-2 ± 0
8	SSP370	2070_2100	4 ± 5.2	7.6 ± 6.9	5.8 ± 12.2	-1.5 ± 9.4	15.9 ± 17.6	-14 ± 16	-2 ± 0
6	SSP585	2030_2060	0 ± 0	4.6 ± 8.1	7.9 ± 11.8	-0.6 ± 10.5	11.9 ± 17.8	-9.9 ± 15.2	-2 ± 0
9	SSP585	2070_2100	2 ± 2.4	8.2 ± 9	11.5 ± 10.8	-1.4 ± 11	20.3 ± 18	-18.3 ± 14.9	-2 ± 0

$EUROPE

EUROPE

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	0.4 ± 1.6	0.7 ± 5.3	0.6 ± 3.5	1.7 ± 6.6	4.1 ± 22.2	-5.8 ± 28.7
2	historical	1970_2000	0 ± 0	0 ± 1.4	0 ± 4.8	0 ± 3	0 ± 5.8	0 ± 22.8	0 ± 28.8
3	historical	1985_2015	0 ± 0	-0.2 ± 1.4	-0.5 ± 4.5	-0.1 ± 2.9	-0.8 ± 5.5	1.1 ± 22.8	-0.3 ± 28.5
4	SSP245	2030_2060	0 ± 0	0.4 ± 1.3	1.4 ± 5	1.2 ± 2.8	3 ± 5.9	7.8 ± 22.3	-10.8 ± 28.4
7	SSP245	2070_2100	0 ± 0	0.7 ± 1.5	1.2 ± 5.9	1.2 ± 3.1	3.1 ± 6.8	10.5 ± 21.6	-13.6 ± 27.6
5	SSP370	2030_2060	0 ± 0	0.6 ± 1.4	0.8 ± 5.8	0.9 ± 3.3	2.3 ± 6.8	7.8 ± 22.5	-10.2 ± 28.6
8	SSP370	2070_2100	0 ± 0	1.6 ± 1.5	2.5 ± 6.8	2.5 ± 4.3	6.6 ± 8.2	7.9 ± 21.5	-14.6 ± 27.6
6	SSP585	2030_2060	0 ± 0	0.8 ± 1.5	1.6 ± 5	1 ± 3.2	3.4 ± 6.1	10.3 ± 22.6	-13.7 ± 28.5
9	SSP585	2070_2100	0 ± 0	1.5 ± 2	2.2 ± 6.4	2.2 ± 3.6	5.9 ± 7.6	12.3 ± 20.6	-18.1 ± 26.7

$EUROPE-AFRICA

EUROPE-AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	-1.1 ± 7.4	-2.9 ± 10.3	-2 ± 13.5	-1.1 ± 4.6	-7.1 ± 19.1	-2.6 ± 21.8	10 ± 7.5
2	historical	1970_2000	0 ± 7.8	0 ± 10.6	0 ± 14.6	0 ± 4.5	0 ± 20.1	0 ± 23.3	0 ± 10.6
3	historical	1985_2015	0.1 ± 8.1	-0.2 ± 10.3	-1.1 ± 12.9	0.4 ± 5.3	-0.8 ± 19.1	0.3 ± 21.8	0.7 ± 10.3
4	SSP245	2030_2060	2.9 ± 8.6	2.4 ± 11.4	3.4 ± 13.6	-0.5 ± 4.7	8.2 ± 20.3	-4.4 ± 23.1	-3.6 ± 7.5
7	SSP245	2070_2100	4.7 ± 7.9	2.1 ± 11.6	3.9 ± 13.1	-0.4 ± 4.5	10.3 ± 19.7	-5.6 ± 21.6	-4.5 ± 7.5
5	SSP370	2030_2060	3.3 ± 8.4	1.7 ± 11.3	2.7 ± 13.5	-0.9 ± 4.3	6.8 ± 19.9	-2.7 ± 21.6	-4 ± 7.5
8	SSP370	2070_2100	6.2 ± 8.3	0.8 ± 11.4	5 ± 14.1	-1.2 ± 4.8	10.8 ± 20.5	-6.3 ± 23.1	-4.5 ± 7.5
6	SSP585	2030_2060	4.1 ± 8.4	2.3 ± 11.2	4.2 ± 13.5	-1.2 ± 4.4	9.4 ± 20	-4.9 ± 21.7	-4.3 ± 7.5
9	SSP585	2070_2100	8.2 ± 7.9	2.6 ± 12.6	4.6 ± 13.5	-1.9 ± 4.3	13.5 ± 20.6	-8.9 ± 21.4	-4.5 ± 7.5

$NORTH-AMERICA

NORTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	0 ± 0.3	1.5 ± 7.1	0.1 ± 3.5	1.6 ± 8	0.6 ± 13.3	-2.3 ± 16
2	historical	1970_2000	0 ± 0	0 ± 0.3	0 ± 6.2	0 ± 3.4	0 ± 7.1	0 ± 13.9	0 ± 16.5
3	historical	1985_2015	0 ± 0	-0.1 ± 0.3	0.3 ± 5.8	-0.1 ± 3.3	0.1 ± 6.7	-0.9 ± 14.2	0.7 ± 16.4
4	SSP245	2030_2060	0 ± 0	0 ± 0.4	2.2 ± 7	0.5 ± 3.5	2.7 ± 7.9	3 ± 13.1	-5.7 ± 16.1
7	SSP245	2070_2100	0 ± 0	0.1 ± 0.4	2.1 ± 6.7	0.9 ± 3.8	3.1 ± 7.7	5.2 ± 13.3	-8.4 ± 16.1
5	SSP370	2030_2060	0 ± 0	0.1 ± 0.4	1.6 ± 6.5	0.4 ± 3.6	2.1 ± 7.5	2.5 ± 13.8	-4.7 ± 16.9
8	SSP370	2070_2100	0.1 ± 0	0.4 ± 0.7	3.3 ± 5.9	1 ± 3.1	4.8 ± 6.8	3.9 ± 14.4	-8.8 ± 16.2
6	SSP585	2030_2060	0 ± 0	0.1 ± 0.4	2.3 ± 6.8	1 ± 3.7	3.4 ± 7.8	3.7 ± 13.4	-7.2 ± 16.5
9	SSP585	2070_2100	0.2 ± 0	0.6 ± 1.4	2.7 ± 7.5	1 ± 3.5	4.5 ± 8.4	7.6 ± 13.5	-12 ± 16.3

$OCEANIA

OCEANIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0.1 ± 0	0.7 ± 37.1	-0.1 ± 24.4	-0.2 ± 4.4	0.5 ± 44.6	-0.4 ± 17.9	-0.2 ± 0.4
2	historical	1970_2000	0 ± 0	0 ± 37.3	0 ± 24.6	0 ± 4.4	0 ± 44.9	0 ± 18.1	0 ± 0.6
3	historical	1985_2015	0 ± 0	2.1 ± 37.1	-1.3 ± 25	0 ± 4.5	0.8 ± 45	-0.9 ± 17.8	0.1 ± 0.5
4	SSP245	2030_2060	5.2 ± 0	4.5 ± 35.3	-5.1 ± 23	-1.1 ± 4	3.5 ± 42.4	-3.4 ± 17.5	-0.2 ± 0.5
7	SSP245	2070_2100	5.3 ± 0	6.6 ± 34.5	-7.7 ± 21.8	-0.6 ± 4.5	3.6 ± 41	-3.1 ± 17.5	-0.4 ± 0.4
5	SSP370	2030_2060	5.5 ± 0	3.5 ± 34.8	-3.4 ± 23.2	-1 ± 4	4.6 ± 42	-4.4 ± 15.1	-0.2 ± 0.4
8	SSP370	2070_2100	5.3 ± 10.2	6.8 ± 32.9	-8.3 ± 22.2	-1.1 ± 4.5	2.7 ± 41.3	-2.5 ± 20.7	-0.2 ± 0.5
6	SSP585	2030_2060	0.4 ± 3.6	4.3 ± 36.6	-3 ± 24.1	-0.1 ± 4.8	1.6 ± 44.2	-1.6 ± 18.8	-0.2 ± 0.5
9	SSP585	2070_2100	4.1 ± 6.4	8.1 ± 34.2	-7.4 ± 21.7	-1.3 ± 4	3.5 ± 41.2	-3.1 ± 17.9	-0.5 ± 0.4

$SOUTH-AMERICA

SOUTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0.1	0.3 ± 2.1	1.3 ± 9.7	0.1 ± 3.5	1.7 ± 10.5	-1.3 ± 13.8	-0.4 ± 2.8
2	historical	1970_2000	0 ± 0.1	0 ± 1.8	0 ± 6.1	0 ± 3.3	0 ± 7.1	0 ± 10.3	0 ± 3.7
3	historical	1985_2015	0.1 ± 0.1	0.2 ± 1.8	0.6 ± 6.7	0.2 ± 3.3	1.1 ± 7.7	-0.9 ± 10.6	-0.2 ± 3.7
4	SSP245	2030_2060	0.1 ± 0.2	0.9 ± 2.6	3.5 ± 10.5	1.5 ± 4.3	6 ± 11.6	-5.1 ± 15	-1 ± 2.8
7	SSP245	2070_2100	0.1 ± 0.1	1.5 ± 3.7	4.4 ± 10.1	2 ± 5.1	8 ± 11.9	-7 ± 16.4	-1 ± 2.7
5	SSP370	2030_2060	0.1 ± 0.1	0.4 ± 2.1	2.4 ± 7	1.5 ± 4.4	4.4 ± 8.5	-4 ± 11.5	-0.4 ± 3.1
8	SSP370	2070_2100	0.1 ± 0.1	1.5 ± 2.8	4.7 ± 8.2	1.4 ± 4.4	7.7 ± 9.8	-7.2 ± 13	-0.5 ± 3.2
6	SSP585	2030_2060	0.1 ± 0.1	1.3 ± 2.8	4.7 ± 10	1.9 ± 5	8 ± 11.5	-6.9 ± 15	-1.1 ± 2.7
9	SSP585	2070_2100	0.5 ± 1	2.6 ± 4	5.7 ± 10.4	2.4 ± 4.5	11.2 ± 12.1	-10.1 ± 14.6	-1.2 ± 2.7

9.3.3 Changes in proportion of aridity categories, by continent, calculated by model

cmip6$cat.AI <- factor(cmip6$cat.AI, levels = c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid", "Cold"))

tab.percent.cont <- cmip6 %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, Continent, cat.AI, source) %>%
  summarise(count = n()) %>%
  ungroup() %>% 
  group_by(period, model, Continent, cat.AI) %>% 
  summarise(mmmean = mean(count, na.rm = T), mmsd = sd(count, na.rm = T)) %>%
  mutate(percent = round(mmmean/sum(mmmean)*100, 1), sd.percent = round(mmsd/sum(mmmean)*100,1)) 
write.table(tab.percent.cont, "tab.percent.continent.txt")

tab.cont <- cmip6 %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, Continent, cat.AI, source) %>%
  summarise(count = n()) %>%
  ungroup() 

list_ref <- list()

for(i in unique(tab.cont$source)){
  tab.cont.i <- filter(tab.cont, source == i)
  tab.cont.i.ref <- filter(tab.cont.i, period == "1970_2000")
  
  for(j in unique(tab.cont.i$Continent)){
    tab.cont.i.j <- filter(tab.cont.i, Continent == j)
    tab.cont.i.j.ref <- filter(tab.cont.i.ref, Continent == j)
    
    for(k in unique(tab.cont.i$cat.AI)){
      
      index <- paste(i, j, k, sep = "_")
      tab.cont.i.j.k <- filter(tab.cont.i.j, cat.AI == k)
      
      lines <- ifelse(dim(tab.cont.i.j.k)[1] > 0, dim(tab.cont.i.j.k)[1] - 1, 1)
      
      tab.cont.i.j.k.ref <- subset(tab.cont.i.j.ref, cat.AI == k)
      
      if(dim(tab.cont.i.j.k.ref)[1] == 0){
        names1 <- names(tab.cont.i.j.k.ref)
        line1 <- c("1970_2000","historical", i, j, k, 0)
        tab.cont.i.j.k.ref <- rbind(tab.cont.i.j.k.ref, line1) %>%
          setNames(names1)
        
      }else{
        tab.cont.i.j.k.ref <- tab.cont.i.j.k.ref
      }
      
      tab.cont.i.j.k.ref <- tab.cont.i.j.k.ref %>% rbind(tab.cont.i.j.k.ref[rep(1,lines),])
      tab.cont.i.j.k$diff <- tab.cont.i.j.k$count - as.numeric(tab.cont.i.j.k.ref$count)  
      
      list_ref[[index]] <- tab.cont.i.j.k
    }
  }
}


tab.cont.diff <- bind_rows(list_ref, .id = "column_label") %>%
  group_by(period, model, Continent, cat.AI) %>% 
  summarise(diff.mean = mean(diff, na.rm = T), diff.sd = sd(diff, na.rm = T))
write.table(tab.cont.diff, "tab.cont.diff.txt")

tab.cont.diff <- read.table("tab.cont.diff.txt")

tab_list_percent <- list()
tab_list_sd <- list()

for(i in unique(tab.cont.diff$Continent)){
  tab.i <- subset(tab.cont.diff, Continent == i)
  
  df.percent <- tab.i %>% reshape2::dcast(period + model ~ cat.AI, value.var = "diff.mean")
  
  df.sd <- tab.i %>% reshape2::dcast(period + model ~ cat.AI, value.var = "diff.sd")
  # sd for a sum: square root of sum of squared sd
  
  drycats <- which(names(df.percent) %in% c("Hyperarid","Arid","Semi-arid","Dry subhumid"))
  
  df.percent <- df.percent %>% mutate("Sum drylands" = rowSums(.[drycats], na.rm = T))
  df.sd <- df.sd %>% mutate("Sum drylands" = round(sqrt(rowSums((.[drycats])^2, na.rm = T)),2))
  
  missing <- setdiff(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"), names(df.percent))
  
  df.percent[, missing] <- 0
  df.sd[,missing] <- 0
  
  df.percent <- df.percent %>% select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))
  df.percent[is.na(df.percent)] <- 0
  
  tab_list_percent[[i]] <- df.percent
  
  df.sd <- df.sd %>% select(c("model","period","Hyperarid","Arid","Semi-arid","Dry subhumid","Sum drylands","Humid","Cold","NA"))
  df.sd[is.na(df.sd)] <- 0
  
  tab_list_sd[[i]] <- df.sd
}

k_list <- list()

for(i in names(tab_list_percent)){
  
  df <- tab_list_percent[[i]] %>% select(c("model", "period"))
  
  df$Hyperarid <- paste(tab_list_percent[[i]]$Hyperarid, tab_list_sd[[i]]$Hyperarid, sep = " ± ")
  df$Arid <- paste(tab_list_percent[[i]]$Arid, tab_list_sd[[i]]$Arid, sep = " ± ")
  df$'Semi-Arid' <- paste(tab_list_percent[[i]][,5], tab_list_sd[[i]][,5], sep = " ± ")
  df$'Dry subhumid' <- paste(tab_list_percent[[i]][,6], tab_list_sd[[i]][,6], sep = " ± ")
  df$'Sum drylands' <- paste(tab_list_percent[[i]][,7], tab_list_sd[[i]][,7], sep = " ± ")
  df$Humid <- paste(tab_list_percent[[i]]$Humid, tab_list_sd[[i]]$Humid, sep = " ± ")
  df$Cold <- paste(tab_list_percent[[i]]$Cold, tab_list_sd[[i]]$Cold, sep = " ± ")
  
  k_list[[i]] <- kable(df[order(df$model, decreasing = F),], caption = i) %>% kable_styling(bootstrap_options = "bordered") %>%
    column_spec(8, italic = T, background = colorvec[5]) %>% row_spec(2, bold = T) 
  
}

k_list

$AFRICA

AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	-5.84615384615385 ± 24.0584224816943	3 ± 22.7449628123093	4.76923076923077 ± 33.9562312547438	-0.384615384615385 ± 11.9620125225542	1.53846153846153 ± 48.91	2.84615384615385 ± 38.6972138399441	0 ± 0
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	-10.5384615384615 ± 19.5474439275974	11 ± 21.6024689946929	5.92307692307692 ± 11.9335768489693	-3.53846153846154 ± 10.1457329669126	2.84615384615388 ± 33.08	-2.69230769230769 ± 15.5155009339296	-2 ± 0
4	SSP245	2030_2060	-10.4615384615385 ± 41.0479706859657	24.3076923076923 ± 29.7388203514772	25.8461538461538 ± 40.2778809311971	2.53846153846154 ± 15.4844835486764	42.2307692307691 ± 66.57	-37.8461538461538 ± 53.1551912702265	0 ± 0
7	SSP245	2070_2100	10.0769230769231 ± 66.8324541153243	9 ± 46.1808040928985	42.0769230769231 ± 51.9157354733956	-2.30769230769231 ± 22.6545970882461	58.8461538461539 ± 99.03	-54.4615384615385 ± 61.6814604569523	0 ± 0
5	SSP370	2030_2060	-29.4615384615385 ± 88.4153223755319	-5.23076923076923 ± 60.9455410539522	45.4615384615385 ± 86.1293749586587	12.9230769230769 ± 45.5621581623419	23.6923076923077 ± 145	-21.3846153846154 ± 43.4291347322249	-30 ± 0
8	SSP370	2070_2100	19.4615384615385 ± 66.2578490754308	18.4615384615385 ± 58.2217247320039	31.0769230769231 ± 46.6341461779198	0.230769230769231 ± 27.395723042578	69.2307692307693 ± 103.47	-67.0769230769231 ± 52.3425600228303	-28 ± 0
6	SSP585	2030_2060	-7.84615384615385 ± 50.8524764291215	13.6153846153846 ± 39.557002038279	33.2307692307692 ± 47.566013507815	0 ± 18.4481254693623	39 ± 82.18	-34.6153846153846 ± 52.663299145829	0 ± 0
9	SSP585	2070_2100	45.5833333333333 ± 209.609792594514	9.25 ± 79.7702667323199	16.5 ± 83.830456604658	-3.83333333333333 ± 20.8232096614826	67.5 ± 240.33	-62.75 ± 98.0381419281663	0 ± 0

$ASIA

ASIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	4.30769230769231 ± 28.0070961703893	8.76923076923077 ± 41.4993049061344	11.3076923076923 ± 61.526938023634	15.8461538461538 ± 39.580816384216	40.2307692307692 ± 88.65	90.2307692307692 ± 176.826258347073	-130.461538461538 ± 287.741091777723
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	-5.15384615384615 ± 19.9743425166977	-12.6153846153846 ± 41.9554097853472	-17 ± 29.8161030317511	-1.61538461538462 ± 17.1975698939243	-36.3846153846154 ± 57.83	-33.6153846153846 ± 169.076185225053	70 ± 234.12461069553
4	SSP245	2030_2060	-15.6923076923077 ± 38.0840837607712	24.2307692307692 ± 40.5301407312176	22.2307692307692 ± 88.9280175630398	27 ± 65.5222608482542	57.7692307692307 ± 123.67	192.692307692308 ± 202.871217202517	-250.461538461538 ± 375.760920308072
7	SSP245	2070_2100	-8.46153846153846 ± 31.9103310977688	23.6153846153846 ± 37.8627751702083	30.0769230769231 ± 101.227023350537	22.8461538461538 ± 72.2286256201401	68.076923076923 ± 133.85	294.461538461538 ± 208.6546650108	-362.538461538462 ± 381.791744145203
5	SSP370	2030_2060	-23.8461538461538 ± 61.3716087370566	21.6153846153846 ± 45.7156764898331	0 ± 55.9359753051052	-2.30769230769231 ± 25.9627445113975	-4.53846153846151 ± 98.28	174.230769230769 ± 120.104228239582	-169.692307692308 ± 100.362164696484
8	SSP370	2070_2100	40.0769230769231 ± 64.3343629517714	25.3076923076923 ± 63.8336178610516	32.6923076923077 ± 53.0383267323179	3.07692307692308 ± 29.0214332820347	101.153846153846 ± 108.94	241.153846153846 ± 179.945753934274	-342.307692307692 ± 156.047740886876
6	SSP585	2030_2060	-1.84615384615385 ± 30.813541811586	23.8461538461538 ± 34.0632895109632	27.0769230769231 ± 91.2044055391163	23.8461538461538 ± 54.9087821661923	72.9230769230768 ± 115.94	250 ± 201.292490338479	-322.923076923077 ± 358.702862533523
9	SSP585	2070_2100	55.0833333333333 ± 178.360441156936	59.5 ± 145.979761984636	77.9166666666667 ± 202.755768142737	22.5 ± 61.9259088089806	215 ± 313.16	324.083333333333 ± 262.604596765602	-539.083333333333 ± 416.369691136338

$CENTRAL-AMERICA

CENTRAL-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	-0.307692307692308 ± 4.67946961710789	-0.538461538461538 ± 17.9704172860815	3.30769230769231 ± 10.7964856580326	2.46153846153846 ± 21.48	-2 ± 28.2105181330179	0 ± 0
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	0 ± 0	1.69230769230769 ± 4.02874288466877	2.53846153846154 ± 6.25320430679856	-0.0769230769230769 ± 4.34859246311452	4.15384615384615 ± 8.62	-4 ± 5.03322295684717	-2 ± 0
4	SSP245	2030_2060	0 ± 0	9.84615384615385 ± 11.1343773500973	17.3076923076923 ± 25.7079255463648	-3.15384615384615 ± 10.0070487977738	24 ± 29.75	-23.5384615384615 ± 26.7507488512483	0 ± 0
7	SSP245	2070_2100	0 ± 0	11 ± 10.2875329080073	22.6153846153846 ± 29.2191103604543	-1.69230769230769 ± 19.0367391088242	31.9230769230769 ± 36.36	-31.4615384615385 ± 27.3484167263097	0 ± 0
5	SSP370	2030_2060	0 ± 0	11.7692307692308 ± 12.3568728929413	3.30769230769231 ± 24.015219533262	0 ± 15.6258333111123	15.0769230769231 ± 31.2	-14.6153846153846 ± 37.2615495775884	0 ± 0
8	SSP370	2070_2100	12 ± 15.556349186104	23.2307692307692 ± 15.0839530967728	18.6923076923077 ± 29.4855462200963	-3.69230769230769 ± 13.4001530798757	50.2307692307692 ± 38.97	-39.6153846153846 ± 23.8381293363471	0 ± 0
6	SSP585	2030_2060	0 ± 0	12.8461538461538 ± 14.0466256555691	21.0769230769231 ± 29.104585945808	-2.46153846153846 ± 16.4195786822489	31.4615384615384 ± 36.25	-31 ± 25.5049014897137	0 ± 0
9	SSP585	2070_2100	6 ± 7.07106781186548	25.1666666666667 ± 24.8845820588256	35.5833333333333 ± 28.5767175840541	-6.16666666666667 ± 18.5709124137396	60.5833333333333 ± 42.79	-55.0833333333333 ± 42.4273889525779	0 ± 0

$EUROPE

EUROPE

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	4.66666666666667 ± 11.211465992321	8.16666666666667 ± 13.2653570238308	7.84615384615385 ± 11.3346530151725	20.6794871794872 ± 20.74	50.3076923076923 ± 106.563271201811	-70 ± 115.036950585453
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	0 ± 0	-2.16666666666667 ± 9.81495457622364	-5.83333333333333 ± 24.5683955643436	-1.23076923076923 ± 8.48679215166962	-9.23076923076923 ± 27.78	12.9230769230769 ± 48.5291347859914	-4.30769230769231 ± 73.0506497979137
4	SSP245	2030_2060	0 ± 0	4.66666666666667 ± 7.41415517572279	17.5 ± 14.5695822614352	11.5 ± 9.56793888700459	33.6666666666667 ± 18.94	97.5384615384615 ± 125.997761477877	-128.538461538462 ± 146.339568233507
7	SSP245	2070_2100	0 ± 0	8.58333333333333 ± 11.1555802044076	21.6153846153846 ± 19.9354728291826	14 ± 9.32737905308882	44.1987179487179 ± 24.68	123.076923076923 ± 120.805257569405	-166.615384615385 ± 139.523318517932
5	SSP370	2030_2060	0 ± 0	7.16666666666667 ± 8.48349556313298	17.6153846153846 ± 11.2215452110249	10.3846153846154 ± 7.83728759987005	35.1666666666667 ± 16.1	90.5384615384615 ± 75.058438771195	-125.153846153846 ± 73.5026599902413
8	SSP370	2070_2100	0 ± 0	19.25 ± 13.7452471124319	37.6923076923077 ± 28.0338623078847	30 ± 31.8093277724213	86.9423076923077 ± 44.57	92.2307692307692 ± 110.037988778235	-177.692307692308 ± 82.173986775403
6	SSP585	2030_2060	0 ± 0	9.83333333333333 ± 8.20014781832768	19.0833333333333 ± 14.2028059540957	12.2307692307692 ± 8.41777728138339	41.1474358974358 ± 18.43	125.923076923077 ± 134.35987839782	-164.846153846154 ± 144.399241776545
9	SSP585	2070_2100	0 ± 0	18.2727272727273 ± 18.9530520449394	37.4166666666667 ± 24.7403179715487	26.0833333333333 ± 24.6298353197659	81.7727272727273 ± 39.72	150.666666666667 ± 131.040821768245	-230.916666666667 ± 148.445556798893

$EUROPE-AFRICA

EUROPE-AFRICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	-2.07692307692308 ± 5.1552810861216	-1.15384615384615 ± 9.04476048923864	1.76923076923077 ± 14.8893353677156	-0.384615384615385 ± 8.18065259762795	-1.84615384615385 ± 19.92	1 ± 20.7163381577601	13 ± 0
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	0.538461538461538 ± 5.82545255202524	0.230769230769231 ± 8.23765587767384	-4.15384615384615 ± 17.8505189179762	2 ± 13.1782649338472	-1.38461538461538 ± 24.37	2.46153846153846 ± 16.7163362444017	-6.5 ± 9.19238815542512
4	SSP245	2030_2060	11.2307692307692 ± 8.13586551586908	7 ± 10.0249688278817	11.6153846153846 ± 19.8432963556061	-3.30769230769231 ± 10.1438373359117	26.5384615384615 ± 25.76	-22.2307692307692 ± 35.5109228410871	-54 ± 0
7	SSP245	2070_2100	19.3846153846154 ± 8.20100056434642	6.84615384615385 ± 11.3860013016434	14.9230769230769 ± 23.0704050624082	-2.53846153846154 ± 10.6584503612187	38.6153846153846 ± 29.03	-36.5 ± 35.8570393748183	0 ± 0
5	SSP370	2030_2060	14 ± 6.7700320038633	6.07692307692308 ± 9.92019437361267	10.5384615384615 ± 14.3097599829358	-4.53846153846154 ± 8.04793332286189	26.076923076923 ± 20.34	-23.0833333333333 ± 20.7121322958894	-56 ± 0
8	SSP370	2070_2100	25.3076923076923 ± 16.7101995568805	-0.923076923076923 ± 23.063179668256	17.6923076923077 ± 22.1260955110499	-6.53846153846154 ± 21.33313301188	35.5384615384615 ± 41.9	-30.9230769230769 ± 28.5115810927815	0 ± 0
6	SSP585	2030_2060	16.9230769230769 ± 9.8696634395635	7.84615384615385 ± 10.1722347089693	16.5384615384615 ± 17.1736978381447	-5.84615384615385 ± 13.0054229320321	35.4615384615384 ± 25.79	-33.1666666666667 ± 32.0137281159003	-57 ± 0
9	SSP585	2070_2100	32.4166666666667 ± 12.273092026005	11.3333333333333 ± 19.90583894847	21.1666666666667 ± 22.0735684528448	-8.66666666666667 ± 18.1174618598308	56.25 ± 36.91	-55.5454545454545 ± 27.4494576863137	0 ± 0

$NORTH-AMERICA

NORTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	0 ± 0	0.375 ± 3.06768875307034	31.75 ± 66.2901679258422	4.92307692307692 ± 27.7172074665467	37.0480769230769 ± 71.92	21.6923076923077 ± 91.341287319759	-56.0769230769231 ± 173.993036995958
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	0 ± 0	0.333333333333333 ± 4.18330013267038	-3 ± 32.2038958682161	-2.69230769230769 ± 36.415163085782	-5.35897435897436 ± 48.79	-19.7692307692308 ± 70.0180379689808	25.0769230769231 ± 113.976943237409
4	SSP245	2030_2060	0 ± 0	3.36363636363636 ± 2.87307247638229	60 ± 104.176132903207	14.1538461538462 ± 34.693530025655	77.5174825174826 ± 109.84	82.1538461538462 ± 94.51705856074	-159.153846153846 ± 214.335735919859
7	SSP245	2070_2100	0 ± 0	4.7 ± 3.05686840482943	56.9230769230769 ± 85.9946331062405	25.4615384615385 ± 49.6548342470556	87.0846153846154 ± 99.35	145.923076923077 ± 122.363162715515	-231.923076923077 ± 215.561847868333
5	SSP370	2030_2060	0 ± 0	5.45454545454545 ± 6.89004552036685	44.4615384615385 ± 45.1582686865787	12.0769230769231 ± 33.2452040111591	61.9930069930071 ± 56.5	70.4615384615385 ± 48.5397009066039	-131.615384615385 ± 77.0427786595846
8	SSP370	2070_2100	4 ± 0	15.6923076923077 ± 20.0307455984736	91.5384615384615 ± 59.6093049009065	28.3846153846154 ± 49.8339550600098	139.615384615385 ± 80.24	108 ± 110.324521299664	-243.923076923077 ± 105.99564577414
6	SSP585	2030_2060	0 ± 0	7.09090909090909 ± 7.86707754448303	62.0769230769231 ± 90.9537075829068	29.8461538461538 ± 46.234630155772	99.013986013986 ± 102.33	103.384615384615 ± 97.9187915754159	-201.307692307692 ± 203.199157074771
9	SSP585	2070_2100	5 ± 0	20.1818181818182 ± 38.4468937154048	91.1666666666667 ± 135.959107042812	30.75 ± 35.4481054757258	147.098484848485 ± 145.67	198.75 ± 99.7835156726801	-339.583333333333 ± 220.157158054023

$OCEANIA

OCEANIA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	1 ± 0	6 ± 42.7668095606862	-1.23076923076923 ± 38.0660694191775	-1.23076923076923 ± 6.80874249468831	4.53846153846154 ± 57.66	-3 ± 12.0623380818148	0.166666666666667 ± 1.16904519445001
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	1 ± 0	5.08333333333333 ± 32.8977985815944	-2.23076923076923 ± 21.2922280271223	1.76923076923077 ± 9.42650382480028	5.62179487179487 ± 40.3	-3 ± 7.85281265959316	-0.25 ± 0.957427107756338
4	SSP245	2030_2060	45 ± 0	50.1538461538462 ± 62.0991762611042	-26.4615384615385 ± 39.9408216085903	-5.84615384615385 ± 12.4354744839522	62.8461538461539 ± 74.87	-19.8461538461538 ± 18.3296033501644	-1.5 ± 0.707106781186548
7	SSP245	2070_2100	45 ± 0	66.6153846153846 ± 81.5041292016734	-48.4615384615385 ± 58.5172558376521	-2.46153846153846 ± 13.0103509087661	60.6923076923076 ± 101.18	-17.3846153846154 ± 18.8172724800845	-5 ± 0
5	SSP370	2030_2060	47 ± 0	42.3846153846154 ± 50.2535877020047	-12.3846153846154 ± 33.3105050035232	-5.30769230769231 ± 10.3068311925038	71.6923076923077 ± 61.17	-27 ± 61.1882341631134	-3.33333333333333 ± 3.21455025366432
8	SSP370	2070_2100	70.6666666666667 ± 84.0376899571456	63.1538461538462 ± 110.812789690425	-57 ± 62.8887907341205	-6.84615384615385 ± 13.452671072109	69.974358974359 ± 153.22	-14.1538461538462 ± 57.5135435554765	-5 ± 0
6	SSP585	2030_2060	23.5 ± 0.707106781186548	33.8461538461538 ± 55.8910340958878	-23.2307692307692 ± 43.1164196220609	-0.769230769230769 ± 12.2348807796524	33.3461538461538 ± 71.65	-12 ± 20.0956048262632	-1.5 ± 0.707106781186548
9	SSP585	2070_2100	60.3333333333333 ± 48.6757982300582	69.8333333333333 ± 126.826892282594	-50 ± 91.7645207731576	-9 ± 14.1806526327567	71.1666666666666 ± 164.55	-23.9166666666667 ± 22.0802599347897	0 ± 0

$SOUTH-AMERICA

SOUTH-AMERICA

	model	period	Hyperarid	Arid	Semi-Arid	Dry subhumid	Sum drylands	Humid	Cold
1	historical	1850_1880	-0.142857142857143 ± 0.377964473009227	5.07692307692308 ± 13.901208211648	20.0769230769231 ± 83.3551253557747	1.30769230769231 ± 18.8121619145019	26.3186813186813 ± 86.58	-20.6153846153846 ± 79.5115698725019	-6.81818181818182 ± 35.0679859182651
2	historical	1970_2000	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0	0 ± 0
3	historical	1985_2015	0.5 ± 1.22474487139159	2.69230769230769 ± 5.87912430248779	10.5384615384615 ± 15.6504280272425	3 ± 9.31844049899624	16.7307692307692 ± 19.18	-14.2307692307692 ± 13.1539211273917	-2.54545454545455 ± 5.75088925929958
4	SSP245	2030_2060	1.57142857142857 ± 1.13389341902768	13.1538461538462 ± 26.6328310982459	55.6153846153846 ± 93.1374776531432	23.4615384615385 ± 24.7810928754679	93.8021978021979 ± 100	-78.4615384615385 ± 97.2167812885335	-18.8 ± 42.023273974734
7	SSP245	2070_2100	1.42857142857143 ± 0.975900072948533	22.6153846153846 ± 49.907478500285	68.9230769230769 ± 90.7647706422684	31 ± 41.3481156362254	123.967032967033 ± 111.53	-107 ± 130.065368180773	-22.7777777777778 ± 43.6886077192264
5	SSP370	2030_2060	1.25 ± 0.886405260427918	5.92307692307692 ± 19.7250329043306	38.1538461538462 ± 31.8821950149979	24.0769230769231 ± 36.5637833984704	69.4038461538462 ± 52.38	-59.4615384615385 ± 55.635143846037	-12.8888888888889 ± 18.6711304186734
8	SSP370	2070_2100	2 ± 2.26778683805536	22.9230769230769 ± 23.8168761541809	73.5384615384615 ± 82.0869613932032	22.6923076923077 ± 35.8849657827727	121.153846153846 ± 92.73	-109.461538461538 ± 109.343964461248	-15 ± 15.1492574075431
6	SSP585	2030_2060	1.5 ± 1.0690449676497	19.2307692307692 ± 29.9197644992789	72.8461538461538 ± 84.535836733942	30.1538461538462 ± 38.093844983685	123.730769230769 ± 97.44	-107.538461538462 ± 99.3576162024628	-20.1 ± 41.1675438503036
9	SSP585	2070_2100	7.375 ± 15.259072242908	40.5 ± 51.561966250681	96.9166666666667 ± 122.811354474571	38 ± 40.7274756488451	182.791666666667 ± 140.12	-160.166666666667 ± 133.095204020781	-29.375 ± 47.7850469737732

9.3.4 Changes in proportion of aridity categories, figure

list_plots <- list()

for(i in names(tab_list_percent)){
  df.percent <- tab_list_percent[[i]]
  df.sd <- tab_list_sd[[i]]
  
  df <- merge(df.percent, df.sd, by = c("model", "period"))
  
  g <- ggplot(data = df)+geom_point(aes(x=period,y= get("Sum drylands.x"), col = model), shape = "\u2605", size = 5, position = position_dodge(width = 0.5))+
    geom_errorbar(aes(x = period, ymin = get("Sum drylands.x")-get("Sum drylands.y"), ymax = get("Sum drylands.x")+get("Sum drylands.y"), col = model), position = position_dodge(width = 0.5))+
    scale_color_manual(values = colorvec[c(1,3,9,11)])+
    ylim(0,100)+
    labs(x="", y = "Percent of drylands", title = i)+
    theme_minimal()
  list_plots[[i]] <- g
}

ggarrange(plotlist = list_plots, ncol = 2, nrow = 4, common.legend = T, legend = "bottom")

tab.flow <- cmip6s %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, Continent, cat.AI) %>%
  summarise(count = n()) %>%
  ungroup() %>% group_by(period, model, Continent) %>% mutate(percent = round(count/sum(count)*100, 1), count = NULL)

df245 <- rbind(subset(tab.flow, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow, model == "SSP245")) %>% mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b245 <- ggplot(data = df245, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  facet_grid(rows = vars(Continent), switch = "y")+
  scale_y_continuous(position = "right")+
  labs(title = "SSP 2-4.5", y = "%", x = "")+
  theme_minimal()

df370 <- rbind(subset(tab.flow, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow, model == "SSP370")) %>% mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b370 <- ggplot(data = df370, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  facet_grid(rows = vars(Continent), switch = "y")+
  scale_y_continuous(position = "right")+  labs(title = "SSP 3-7.0", y = "%", x = "")+
  theme_minimal()


df585 <- rbind(subset(tab.flow, model == "historical" & period %in% c("1850_1880", "1970_2000")), subset(tab.flow, model == "SSP585")) %>% mutate(lab.y = (rev(cumsum(rev(percent)))) - percent*0.5) %>% 
  subset(cat.AI %in% c("Hyperarid", "Arid", "Semi-arid", "Dry subhumid", "Humid","Cold"))

b585 <- ggplot(data = df585, aes(x = period, y = percent))+
  geom_col(aes(group = cat.AI, col = cat.AI, fill = cat.AI))+
  geom_label(aes(y = lab.y, label = percent))+
  scale_color_manual(values = col.cat, aesthetics = c("col", "fill"), na.translate = F)+
  facet_grid(rows = vars(Continent), switch = "y")+
  scale_y_continuous(position = "right")+  labs(title = "SSP 5-8.5", y = "%", x = "")+
  theme_minimal()


ggarrange(plotlist = list(b245, b370, b585), common.legend = T, legend = "bottom", ncol = 3)

tab.percent.region <- cmip6 %>% subset(!Continent %in% c("SOUTHERN","PACIFIC","ATLANTIC","INDIAN","ARCTIC") & lat > -55 ) %>%
  group_by(period, model, Continent, Name, AI) %>%
  summarise(AI.mean = mean(AI), AI.sd= sd(AI)) %>%
  ungroup()

  
ggplot(subset(tab.percent.region, Continent == "AFRICA"))+
  geom_point(aes(x=period, y = AI.mean, col = Name), position = "jitter", alpha = 0.5)

10 Changes of variables over time and scenarios

10.1 Temperature

10.1.1 Absolute temperature

quantile(cmip6s$t.mean, probs = seq(0,1,0.1)) 
##          0%         10%         20%         30%         40%         50% 
## -53.4587313 -42.5828261 -28.4644864 -17.7502317  -7.2080554  -0.2580962 
##         60%         70%         80%         90%        100% 
##   7.1198406  15.3318920  22.8072968  26.4971129  34.3727997

t.breaks <- c(-50,-40,-30,-20,-10,0,10,20,30)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#fddbc7", "#f4a582", "#d6604d", "#b2182b", "#67001f")

map_list <- list()

for(i in c("1850_1880", "1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = t.mean))+
    borders(colour = "grey60")+
    binned_scale(aesthetics = "fill", breaks = t.breaks, palette = function(x) colscale,
                 guide = guide_legend(label.theme = element_text(angle = 0)))+
    labs(title = i, fill = "Mean annual temperature, °C")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right")
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = t.mean))+
      borders(colour = "grey60")+
      binned_scale(aesthetics = "fill", breaks = t.breaks, palette = function(x) colscale,
                   guide = guide_legend(label.theme = element_text(angle = 0)))+
      labs(fill = "Mean annual temperature, °C", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}



ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

#labels = c("2030-2060", "2070-2100"

10.1.2 Temperature anomalies

quantile(cmip6s$diff.t.mean, probs = seq(0,1,0.1)) 
##         0%        10%        20%        30%        40%        50%        60% 
## -3.5775507 -0.5616741  0.0000000  1.1670755  1.8544034  2.2566573  2.6600679 
##        70%        80%        90%       100% 
##  3.1444949  3.7454555  4.5875606 10.7627615

t.breaks <- c(-3,-0.5,0,1.5,2,3,4,5,10)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#fddbc7", "#f4a582", "#d6604d", "#b2182b", "#67001f")

map_list <- list()

for(i in c("1850_1880", "1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = diff.t.mean))+
    borders(colour = "grey60")+
   scale_fill_gradient2(low = "#3182bd",mid = "white", high = "#b2182b", midpoint = 0)+
    labs(title = i, fill = "Temperature anomalies compared to 1970-2000, °C")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right", legend.text =  element_text(angle = 0))
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_tile(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = diff.t.mean))+
      borders(colour = "grey60")+
   scale_fill_gradient2(low = "#3182bd",mid = "white", high = "#b2182b", midpoint = 0)+
      labs(fill = "Temperature anomalies compared to 1970-2000, °C", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom", legend.text =  element_text(angle = 0))
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

10.2 Precipitation

10.2.1 Absolute values of precipitations

quantile(cmip6s$pr.mean, probs = seq(0,1,0.1)) 
##          0%         10%         20%         30%         40%         50% 
##    8.579696   55.000319  104.132310  200.852004  337.042710  463.946055 
##         60%         70%         80%         90%        100% 
##  579.741473  712.099280  959.490392 1431.315602 6863.723388

p.breaks <- c(10,50,100,200,300,400,500,600,700)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#fddbc7", "#f4a582", "#d6604d", "#b2182b", "#67001f")

map_list <- list()

for(i in c("1850_1880","1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = pr.mean))+
    borders(colour = "grey60")+
    binned_scale(aesthetics = "fill", breaks = p.breaks, palette = function(x) rev(colscale),
                 guide = guide_legend(label.theme = element_text(angle = 0)))+
    labs(title = i, fill = "Annual mean\nprecipitation, mm")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right")
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = pr.mean))+
      borders(colour = "grey60")+
      binned_scale(aesthetics = "fill", breaks = p.breaks, palette = function(x) rev(colscale),
                   guide = guide_legend(label.theme = element_text(angle = 0)))+
      labs(fill = "Annual mean\nprecipitation, mm", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom")
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

10.2.2 Precipitation anomalies

cmip6s$pr.percent <- with(cmip6s, diff.pr.mean/spr.ref.mean*100)

quantile(cmip6s$pr.percent, probs = seq(0,1,0.1)) 
##           0%          10%          20%          30%          40%          50% 
## -42.02999493  -2.87352237   0.00000000   0.04663897   2.46136429   5.19906939 
##          60%          70%          80%          90%         100% 
##   9.34748067  13.67239987  18.11758968  24.02085190 148.15226168

p.breaks <- c(-43,-3,0,3,5,9,13,18,24,148)
colscale <- c("#08519c", "#3182bd", "#6baed6", "#9ecae1", "#c6dbef", "#fddbc7", "#f4a582", "#d6604d", "#b2182b", "#67001f")

map_list <- list()

for(i in c("1850_1880","1970_2000","1985_2015")){
  g <- ggplot() + geom_raster(data = subset(cmip6s, period == i & model == "historical"), aes(x=lon, y = lat,  fill = pr.percent))+
    borders(colour = "grey60")+
   scale_fill_gradient2(high = "#3182bd",mid = "white", low = "#b2182b", midpoint = 0)+
    labs(title = i, fill = "Precipitation anomalies in %")+
    theme_void()+ylim(-55,90)+
    theme(legend.position = "right", legend.text =  element_text(angle = 0))
  map_list[[i]] <- g
}

ggpubr::ggarrange(plotlist = map_list, ncol = 1, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = pr.percent))+
      borders(colour = "grey60")+
   scale_fill_gradient2(high = "#3182bd", mid = "white", low = "#b2182b", midpoint = 0)+
      labs(fill = "Precipitation anomalies in %", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom", legend.text =  element_text(angle = 0))
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")

map_list <- list()

for(i in c("SSP245","SSP370","SSP585")){
  for (j in c("2030_2060","2070_2100")){
    index <- paste(i, j, sep = " , ")
    g <- ggplot() + geom_raster(data = subset(cmip6s, period == j & model == i), aes(x=lon, y = lat,  fill = pr.percent))+
      borders(colour = "grey60")+
   scale_fill_binned_divergingx(mid = 0, palette = "RdBu", rev = F, n_interp = 15)+
      labs(fill = "Precipitation anomalies in %", title = index)+
      theme_void()+ylim(-55,90)+
      theme(legend.position = "bottom", legend.text =  element_text(angle = 0))
    map_list[[index]] <- g
  }}

ggpubr::ggarrange(plotlist = map_list,ncol = 2, nrow = 3, common.legend = T, legend = "bottom")